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ABSTRACT 

Two  very  important  and  fundamental  topics  in  the  theory 
of  dynamical  systems  are  observability  and  identif iability . 
There  has  been  limited  work  on  the  observability  and  iden- 
tif iability  of  nonlinear  dynamical  systems.  This  problem 
is  addressed  by  viewing  the  observation  process  as  a non- 
linear mapping  of  the  initial  states  (and  parameters)  to  the 
output  measurements.  The  problem  of  observability  and  iden- 
tifiability  is  then  reduced  to  finding  conditions  under 
which  the  nonlinear  mapping  is  one-to-one.  Using  such  an 
approach,  the  theory  of  nonlinear  functions  is  extended  to 
find  new  sufficient  conditions  for  local  and  global  observ- 
ability and  identif iability  for  nonlinear  dynamical  systems. 
The  nonlinear  theory  is  then  applied  to  the  problem  of 
determining  the  identif iability  of  the  optimal  control  model 
for  the  human  operator.  The  model  assumes  that  the  operator 
performance  is  dictated  by  the  desire  to  behave  optimally 
with  respect  to  a chosen  cost  functional  under  constraints. 
Within  the  model  structure  are  a number  of  parameters,  the 
values  of  which  determine  the  model  response.  Previous 
attempts  to  establish  the  identif iability  of  the  model 
parameters  have  been  limited  to  linear  system  theory.  The 
identif iability  of  the  model  parameters  is  established 
using  the  previously  developed  nonlinear  theory,  a gradient 
computational  technique  is  used  to  estimate  model  parameters, 
and  the  results  of  model  response  to  experimental  simulator 
data  are  presented. 
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Chapter  I 
INTRODUCTION 

1.1  Background 

In  this  research,  we  will  be  concerned  with  certain 
properties  of  dynamical  systems.  By  a dynamical  system,  we 
mean  a structure  which  at  each  moment  of  time  t e [tQ,  tf] 
receives  some  input  u(t)  and  emits  an  observable  output 
y(t).  It  is  noted  that  the  important  class  of  models 
referred  to  as  dynamical  systems  are  more  formally  defined 


in  the  literature  (e.g.,  see  Ref  17).  Knowledge  of  the  struc- 
ture of  a dynamical  system  allows  one  to  formulate  a mathe- 
matical model  of  the  system  which  includes  a state-transition 
function  (or  solution  or  trajectory)  g(t,  x^,  u,  4>)  whose 
value  is  the  state  x(t)  resulting  at  time  t from  the  initial 

state  Xq#  system  parameter  vector  ((),  and  input  function  u. 

It  is  noted  that  the  state  at  time  t is  dependent  not  only  on 
the  current  value  of  u,  but  also  the  past  history  of  u from 
t^  to  t.  The  observed  output  or  observation  process  is 


described  by  y(t)  = h(t,  x(t),  <j))  for  all  t e tf] 

continuous  observations  or  y(t^^)  = h(tj^,  x(t^),  <{))  with 
tji^  e C^o'  ^f^  discrete  observations.  Two  very  important 
and  fundamental  topics  in  the  theory  of  dynamical  systems 
are  observability  and  identifiability.  A dynamical  system 
is  observable  if,  with  knowledge  of  y(t)  and  u(t)  on 
[tof  tjl  , one  can  uniquely  determine  the  initial  state,  x^. 
If  Xq  is  found,  then  the  state  vector  x(t)  can  be  computed 
using  the  state-transition  function  for  any  time  t e 


i 
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Observability  is  an  important  property  since  control  of 
dynamical  systems  is  ordinarily  implemented  on  the  basis  of 
knowledge  of  the  system  state  x(t).  The  related  problem  of 
system  identif lability  is  concerned  with  determining  a 
unique  parameter  vector  (p  with  knowledge  of  y(t)  and  u(t) 
on  Ct^,  t^].  Identif lability  is  a fundamental  property 
since  one  wants  to  select  a model  structure  and  an  input/ 
output  relation  which  will  allow  the  unique  determination 
of  the  unknown  parameter  values. 

Conditions  for  observability  of  linear  systems  have 
been  studied  extensively.  Some  original  efforts  concerning 
observability  were  due  to  Kalman  (Refs  18,  19) , Kreindler 
and  Sarachik  (Ref  30),  and  Gilbert  (Ref  10).  Likewise  many 

researchers  have  investigated  identification  of  linear 

' • 

systems.  For  example.  Bellman  and  Astrom  (Ref  6),  Mehra 
(Ref  36),  Denham  (Ref  8),  and  others  investigate  conditions 
under  which  certain  canonical  forms  of  linear  time  invar- 
iant systems  will  be  identifiable.  Glover  and  Williams 
(Ref  11)  develop  conditions  for  both  local  and  global 
identifiability  from  the  system  frequency  domain  transfer 
function.  Tse  and  Anton  (Ref  50),  among  others,  have  formu- 
lated necessary  and  sufficient  conditions  for  assessing 
identifiability  for  linear  systems  with  stochastic  disturb- 
ances and  observation  noise. 

The  investigation  of  observability  and  identifiability 
of  nonlinear  systems  has  been  limited.  Kostyukouskii 
(Refs  28,  29)  and  Griffith  and  Kumar  (Ref  17)  report 
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observability  conditions  for  nonlinear  systems  under  con- 
tinuous observation  processes.  Grewall  and  Glover  (Ref  12) 
address  local  identifiability  by  linearizing  the  state 
equations  in  the  neighborhood  of  interest. 

The  approach  taken  in  this  research  is  to  view  the 
observation  process  as  a nonlinear  mapping  of  the  initial 
states  and  parameters  to  the  observed  output.  For  example, 
for  a discrete  observation  process  the  initial  states 
are  mapped  via  the  state  transition  and  output  functions 
to  the  observation  sequence,  which  may  be  denoted  by  an 
output  vector  Y.  That  is  Y(Xq)  may  be  viewed  as  a non- 
linear mapping,  Y : D C -*•  r"*.  The  observability  problem 
is  then  reduced  to  finding  conditions  v/here  the  mapping  Y 
is  one-to-one.  The  identification  problem  can  be  approached 
in  a similar  manner.  With  this  approach,  the  theory  of  non- 
linear functions  can  be  extended  to  find  sufficient  con- 
ditions for  local  and  global  observability  and  identifi- 
ability for  nonlinear  dynamical  systems. 

A practical  example  of  a case  where  a model  structure 
for  dynamical  systems  i^  defined,  but  where  the  question  of 
identifiability  of  model  parameters  is  unresolved,  is  the 
optimal  control  model  for  the  human  operator  (Refs  20,  21, 
22,  23,  25).  The  optimal  control  model  assumes  that  oper- 
ator performance  is  dictated  by  the  desire  to  behave  opti- 
mally with  respect  to  a chosen  cost  functional.  Using  the 
techniques  of  modern  control  and  estimation  theory  and 
incorporating  factors  to  account  for  inherent  human 


limitations  results  in  a model  structure  shovm  in  Fig.  1.1. 
Within  this  model  structure  are  a nvunber  of  parameters,  the 
values  of  which  determine  the  model  response. 
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Given  the  basic  structure  of  the  optimal  control  model, 
one  needs  to  address  the  question  of  identifiability  of  the 

1 

model  parameters  and,  if  identifiable,  use  experimental  data 
to  estimate  their  values.  Previous  attempts  to  establish 
the  identifiability  of  the  model  parameters  have  been 
limited  to  linear  system  theory  (Refs  40,  41).  In  this 
research  the  identifiability  of  the  model  will  be  analyzed 
using  the  theory  for  inear  dynamical  systems  as  previ- 

ously discussed.  The  ncr'  near  theory  is  required  since, 
without  simplifying  assumptions,  the  model  input/output 
mathematical  description  is  nonlinear  with  respect  to  the 
unknown  parameters.  The  model  is  used  in  a tracking 
function  to  correspond  to  an  available  simulator  system  and 
the  output  (or  measured  data)  consists  of  ensemble  averages 
of  many  time  histories  of  tracking  errors.  The  model  mean 
and  covariance  predictions  will  be  developed  in  a manner 
similar  to  Kleinman's  work  on  modeling  a tracking  problem 
(Refs  20,  23).  The  use  of  ensemble  averaged  data  as  the 
measured  output  for  establishing  identifiability  differs 
from  the  more  usual  approach  of  using  the  time  histories 
from  individual  trials  (Refs  44,  48,  49).  However,  assuming 
this  type  of  experimental  data  is  available  and  the  system 
is  identifiable  under  this  type  of  observation  process,  it 
is  appealing  to  adjust  system  parameters  to  match  the 


r 


average  response  over  many  trials  and  different  operators. 
After  the  identif iability  of  the  optimal  control  model  is 
established  for  a specified  input/output  arrangement,  a 
computational  procedure  will  be  used  to  adjust  the  model 
parameters  in  an  attempt  to  minimize  the  mean  square  differ- 
ence between  experimental  output  data  and  model  output  data. 
This  differs  from  previous  work  where  nominal  values  have 
been  established  for  the  model  parameters  using  a subjective 
criteria  to  match  experimental  to  model  data  and  a heuristic 
procedure  for  adjusting  parameters  (Refs  20,  23). 

1.2  Objectives  and  Organization 

The  first  objective  of  this  research  is  to  develop 
sufficient  conditions  for  establishing  observability  and 

I 

E 

I identif iability  of  nonlinear  dynamical  systems.  In 

I Chapter  II  the  observability  question  is  addressed  and 

1 

sufficient  conditions  are  developed  for  local  and  global 
observability  of  nonlinear  systems.  The  identification 
problem  is  analyzed  as  a special  case  of  observability  in 
Chapter  III,  so  that  the  results  of  this  chapter  include 
, both  local  and  global  conditions  for  identifiability  of  non- 

I linear  systems.  In  both  Chapter  II  and  Chapter  III, 

I"  examples  are  presented  to  illustrate  the  application  of  the 

' theory  to  nonlinear  systems,  as  well  as  application  to 

linear  systems.  In  Chapter  IV,  some  theorems  concerning 

j 

1 minimization  of  cost  functionals  are  developed,  and  the 


t 

relation  of  the  problem  of  minimizing  a least  squares  cost 

6 


functional  to  identifiability  is  presented.  In  addition, 
conqputational  procedures  which  can  be  used  to  estimate 
model  parameters  are  discussed. 

A second  major  objective  is  to  use  the  theory  developed 
for  assessing  nonlinear  dynamical  systems  to  address  the 
identifiability  of  the  optimal  control  model  for  human  oper- 
ators. Chapter  V develops  the  mathematical  description  of 
the  model,  including  the  state  space  equations  for  propa- 
gating mean  states  and  covariances.  Then  in  Chapter  VI, 
the  simulator,  from  which  experimental  data  is  obtained,  is 
described  and  the  mathematical  description  developed  in 
Chapter  V is  placed  in  the  context  of  the  simulator.  In 
Chapter  VII,  the  identifiability  of  the  model  parameters  is 
established  using  the  nonlinear  theory  and  the  input/output 
relationships  developed  in  Chapters  V and  VI.  Also  the 
results  of  using  a gradient  computational  technique  to  esti- 
mate model  parameters  are  presented. 

A third  objective  is  to  provide  further  insight  into 
the  modeling  of  human  response.  The  research  is  unique 
with  respect  to  system  ^namics,  sensor  displays,  and  con- 
trol mechanisms.  Other  aspects,  such  as  operator  thresh- 
olds, and  other  possible  model  refinements,  are  discussed 
in  Chapters  V,  VI,  and  VII  as  the  particular  aspect  of  the 
model  arises. 

I 

Finally,  Chapter  VIII  summarizes  the  results  and  their 
significance,  and  presents  recommendations  for  further 
research. 
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Chapter  II 


OBSERVABILITY  OF  NONLINEAR  DYNAMICAL  SYSTEMS 

The  observability  problem  is  concerned  with  determining 
conditions  under  which  knowledge  of  the  input  and  observed 
output  data  uniquely  determine  the  state  of  the  system. 

This  problem  has  been  discussed  extensively  for  linear 
systems;  however,  the  investigation  of  the  nonlinear  problem 
is  limited.  Kostyukovskii  (Refs  28,  29)  and  Griffith  and 
Kumar  (Ref  13)  report  observability  conditions  for  nonlinear 
systems  under  continuous  observation  processes  on  a time 
interval  Grewal  and  Glover  (Ref  12)  address  the 

related  problem  of  local  identifiability  of  nonlinear 
systems  by  linearizing  the  state  equations  in  the  neighbor- 
hood of  interest.  The  approach  taken  in  this  chapter  differs 
from  any  of  these  investigations.  Totally  new  results  are 
obtained  for  determining  the  observability  of  nonlinear 
dynamical  systems.  Results  are  derived  first  for  discrete 
observations  and  then  extended  to  continuous  observations. 
Examples  are  presented  and  the  relation  of  the  theory  to 
linear  systems  is  illustrated.  In  Chapter  III,  this  theory 
will  be  applied  to  parameter  identification  in  nonlinear 
dynamical  systems. 
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2.1  Nonlinear  State  Equations  and  Definition  of 
observability 

We  consider  nonlinear  systems  described  by 

x(t)  = f(t,  x(t),  u(t))  " *o  (2.1.1) 

for  all  t e CtQ,  tf]  and  the  discrete  observation  process 

y(t^)  = h(t^,  x(tj^))  i = 1,  2,  • • • k (2.1.2) 


where  x(t)  is  an  n vector,  y(t)  is  an  m vector,  the  input  u 
is  an  r vector  valued  function  of  t,  and  t.  eft  , t J . 
Assume  f (•,»,• ) has  continuous  partial  derivatives  with 
respect  to  its  first  two  arguments,  u is  a member  of  a set 
U of  admissible  input  functions,  where  U C C the  space  of 
continuous  functions  on  ^f^  * is  continuous  and 

has  continuous  partial  derivatives  with  respect  to  both  of 
its  arguments. 

A solution  to  Eg  (2.1.1)  is  denoted  by  g(t,  x^,  u) 
for  all  t e Et^,  t^].  For  an  admissible  input  u,  we  are 
interested  in  finding  sufficient  conditions  to  guarantee  a 
one-to-one  correspondence  between  the  initial  conditions  x^ 
and  the  observed  output*sequence  (y(tj^)}.  y(t^)  is  an  m 

vector  denoting  the  discrete  observations 


y(tj^)  = h(tj^,  g(tj^,  x^,  u) ) i = 1,  2,  • • • k (2.1.3) 


This  observation  process  can  be  viewed  as  one  large  obser- 
vation vector  of  dimension  mk.  Let  Y denote  this  vector  so 


y = Cy'"(ti)  y"'(t2)  • • • y'^(tk)]' 


(2.1.4) 


i ^ ' 


that 


The  nonlinear  system  is  observable  under  the  observation 


process  Y with  the  input  u if  there  is  a one-to-one  corre- 
spondence between  the  initial  conditions  and  the  obser- 
vation vector  Y(x^).  V7e  say  the  system  is  locally  observ- 
able at  Xq  when  the  one-to-one  correspondence  holds  in  a 
neighborhood  of  Xq*  When  there  is  a one-to-one  corre- 
spondence between  x^  and  Y for  all  x^  £ 

system  is  globally  observable  on  under  the  observation 
process  Y with  input  u. 

2.2  Local  Observability  of  Nonlinear  Dynamical  Systems 

In  this  section,  we  will  develop  theorems  which  lead 
to  sufficient  conditions  for  a nonlinear  system  to  be 
locally  observable. 


Theorem  2 . 1 

Consider  the  nonlinear  function  y = g(x)  where  y is  a 

k vector  and  x is  an  n vector  with  k ^ n.  This  is  denoted 
n k 

by  g : X CR  -»■  R . The  Jacobian  matrix  will  be  denoted  by 
g'(x)  and  is  defined  as  follows: 


Let  be  a point  in  X which  is  mapped  by  g to  and  where 
g(x)  has  a strong  Frechet  derivative  at  x^  and  the  k x n 
Jacobian  matrix  g'(x^)  has  rank  n.  If  we  denote  g(x^)  by  y, 
then  .an  inverse  function  g ^ exists  which  maps  points  in  a 
neighborhood  of  y^  to  a neighborhood  of  x^;  that  is  g is 
locally  one-to-one  at  x^. 


Proof:  For  k = n,  this  is  a statement  of  the  inverse 

function  theorem  (see,  e.g..  Ref  39:125).  For  this  special 
case  of  k = n with  g : X Cr”  R*^,  g is  termed  a local 
homeomorphism  at  x^.  Thus,  at  x^  there  is  an  open  neighbor- 
hood 6 e X such  that  when  g is  restricted  to  6,  the  mapping 
is  one-to-one  and  continuous. 

For  k > n,  if  g' (x)  has  a rank  of  n,  there  is  at  least 
one  n X n submatrix  of  g' (x^)  (say  g^ (x^) ) such  that 
lgi(Xo)l  f 0.  Thus  we  can  define  a locally  homeomorphic 
mapping  g^^  : 6 CR^  -*■  R*^  where  6 is  an  open  neighborhood  of 
Xq.  Further  let  g^  (Xq)  equal  the  n components  of  g(x^) 
corresponding  to  the  rows  of  the  submatrix  gj^  (x^)  which  is 
nonsingular.  If  there  is  more  than  one  n x n submatrix  of 
g' (Xq)  which  is  nonsingular,  arbitrarily  choose  the  non- 
singular matrix  incorporating  the  uppermost  rows  so  that  a 
unique  g^^  is  obtained.  Now  suppose  g restricted  to  6 were 
not  one-to-one.  Then  there  exists  Xj^  e 6 and  X2  e 6 with 

X2  and  g(Xj^)  = g(X2).  But  this  contradicts  g^^  being  one- 
to-cr.'  at  Xq  since  by  the  definition  of  g^^,  g(xj^)  = g(x2) 
implies  gj^(Xj^)  = gj^(x2). 

Q.E.D. 


Theorem  2.1  requires  that  g(x)  have  a strong  Frechet  deriv- 


ative at  Xq.  Appendix  A defines  a Frechet  derivative,  a 
strong  Frechet  derivative,  and  a Gateaux-derivative,  as 
well  as  some  useful  properties  of  these  derivatives. 
Throughout  this  report,  a derivative  means  in  the  Frechet 
sense  except  where  the  phrase  "continuously  differentiable" 
is  used.  If  g is  continuously  differentiable  on  an  open 
set  Dq,  this  means  g has  a continuous  Gateaux-derivative  on 
Dq.  From  Appendix  A,  we  note  that  a continuous  Gateaux- 
derivative  at  Xq  implies  a continuous  Frechet  derivative. 

We  also  note  that  if  g has  a Frechet  derivative  at  each 
point  of  an  open  neighborhood  of  x^,  then  g*  is  strong  at 
Xq  if  and  only  if  g'  is  continuous  at  x^.  Thus  in  all 
sufficiency  theorems,  we  could  replace  "strong  derivative 
at  Xq"  by  "has  a derivative  in  a neighborhood  of  x^  which 
is  continuous  at  x^."  However,  we  should  note  that  g can 
have  a strong  derivative  at  x^^  even  though  g is  not  differ- 
entiable at  all  points  of  an  open  neighborhood  of  Xq. 

Theorem  2 . 2 

Consider  the  nonlinear  system  described  by  Eq  (2.1.1) 
with  input  u and  initial  condition  x^.  The  discrete  obser- 
vation process  Y is  defined  by  Eq  (2.1.4).  Assume  the 
dimension  of  Y is  equal  to  or  greater  than  the  dimension  of 
x(t) ; i.e.,  mk  ^ n.  If  Y(x)  has  a strong  derivative  at  Xq 
and  the  rank  of  the  mk  x n Jacobian  matrix  3Y(xq)/3xq  is 
equal  to  n,  then  the  nonlinear  system  is  locally  observable 
at  Xq  under  the  observation  process  Y with  input  u. 
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Proof:  From  Theorem  2.1,  an  inverse  function  exists  which 

maps  points  in  a neighborhood  of  Y(x^)  to  a neighborhood  of 
Xq.  Thus  Y is  locally  one-to-one  at  x^. 

Q.E . D. 

Although  Theorem  2.2  gives  sufficient  conditions  for 
local  observability,  the  following  lemma  leads  to  conditions 
which  may  be  more  practical  to  apply. 

Lemma  2 . 1 

A necessary  and  sufficient  condition  that  vectors 
*1'  *2,  • • • , Xjj  be  linearly  independent  is  that  their 
Gram  matrix  be  nonsingular.  The  Gram  matrix  G is  defined 
in  terms  of  the  inner  products  of  Xj^  and  x j , (x^^,  Xj)  as 
follows ; 


<=‘l'  >‘2>  ■ ■ ■ 

(Xj^, 

G = [ (x^  , Xj  ) ] = 

(X2, 

• 

• 

xi) 

(X2,  X2)  • • • 

(X2f 

• 

• 

x„) 

• 

• 

j^n' 

^1^ 

(Xn,  Xj)  ■ • • 

'*n' 

= Cx,x^ 
^12 

• • 

• [’‘1*2  • • • 

Proof;  See  Ref  38:378. 

Theorem  2.3 

Let  Y (x^)  have  a strong  derivative  at  x^.  A sufficient 
condition  for  a nonlinear  system  to  be  locally  observable 
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at  X under  the  observation  process  Y with  input  u is  that 
o 


k |3y(t.)  " r 9y(t. ) 

y i_  i_ 

i=l  ^^o  9^0 


be  nonsingular. 


Proof;  From  Theorem  2.2,  we  know  a sufficient  condition  is 
that  the  rank  of  the  mk  x n matrix  be  equal  to  n.  Thus  it 


is  sufficient  that  the  columns  of  3Y(x^)/9x^  be  independent. 


Recall  that  Y = Cy  (t-j^)y  (t2) 


(t,  )]  , so  that 


9y(tj^)  9y(tj^) 


3y  (t]^)' 


9Y(x^) 


9y (t2) 
9x  . 


3y (t2) 
9x  _ 


9y  (t2> 


9y(t  ) 
k 


3y(t  ) 
k 


9y(t  ) 
k 


mk  X n 


(2.1.5a) 


V 


1 


9Xn 


pyi(tj)- 

on 

• 

• 

• 

• 

3ym<^l) 

9x^, 

9x 

ol 

• 

on 

• 

• 

• 

• 

• 

pYl(t^)1 

■3yi(t^)-| 

^^ol 

• 

• 

^*Ol 

• 

• 

• 

m k 

• • • 

• 

m K 

3^01 

on 

(2.1.5b) 


From  Lemma  2.1,  a necessary  and  sufficient  condition  that 
the  columns  of  9Y(Xq)/9x^  be  independent  is  that  the  n x n 
Gram  matrix  of  9Y(Xq)/9Xq  (denoted  by  ^ non- 

singular. But  recall  from  the  definition  of  the  Gram  matrix 
that 


'9Y/9x, 


-.T  ^ 


r 

o 

■ 9Y(x^)‘ 

9Xq  -- 

(2.1.6) 


Substituting  Eq  (2.1.5)  into  Eq  (2.1.6)  yields 


k 


9y (t  , X ) 
1 o 


,-iT  r 


9x 


9y (t . , X ) 
1 o 


3x^ 


(2.1.7) 


Q . E . D . 


L 


Note  that  in  Theorems  2.2  and  2.3  we  have  placed  no  require- 
ments on  the  uniqueness  of  the  solution  of  the  system  Eq 


I 


(2.1.1).  Of  course  the  solution  as  mapped  by  Y (x^)  must  be 
unique  in  a local  neighborhood  of  x^,  but  this  is  assured 
by  the  conditions  of  Theorems  2.2  and  2.3. 

At  this  point,  an  example  will  be  presented  to  illus- 
trate the  use  of  the  previous  theorems.  Let 


Xj^(t)  = Xj^(t) 


X2(t)  = 2x2 (t) 


with 


xi(0)  = x^^ 


X2(0)  = x^2 


for  t e [0,t^].  Thus 


= "'ol® 


The  observations  are 


y(t.)  = x^{t^)  + X ^t.) 


i = 1,  2,  3 


so  that  the  observation.^rocess  is  described  by 


Y = 


(t^)  + 

* =<2^  <'=2* 
* >=2^(t3) 


16 
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Therefore,  we  can  state  that  3Y(x^)/3x^  has  a rank  of 


two  for  Xq  e {R^  ~ ^*o1  “ ® ^o2  ” where  {A  ~ B} 

denotes  the  set  of  all  points  in  A which  are  not  in  B. 
Prpm  Theorem  2.2,  we  can  conclude  that  the  system  is 
locally  observable  under  the  observation  process  Y for 
Xj,  e {R^  ~ {xqi  = 0 or  x^2  = 


If  we  were  to  use  Theorem  2.3,  we  first  form 


0 or 


O 

Again  this  is  nonsingular  for  e {R  - = 

Xq2  “ 0}},  so  that  from  Theorem  2.3,  the  system  is  locally 
observable  on  this  set  under  the  observation  process  Y. 

Theorems  2.2  and  2.3  can  be  used  to  make  statements  of 
sufficient  conditions  for  observability  of  nonlinear  systems 
under  continuous  observations. 

Theorem  2 . 4 

Let  y(t)  = h(t,  x(t))  be  the  continuous  observation 
process  associated  with  the  nonlinear  system 

x(t)  = f(t,  x(t),  u(t))  x(t  ) = X 

o o 

for  all  t e • 

Let  Y be  the  discrete  observation  process  vector  formed  by 

n 

the  components  y(tj^),  i = l,  2***k  with  k ^ m where 
ti  e [t^f  tjl  are  any  k time  points  on  t^l . If  Y has 

a strong  derivative  at  x^  and  the  rank  of  the  Jacobian 
9Y(Xq)/3Xq  formed  by  any  such  k time  points  is  equal  to  n, 
then  the  nonlinear  system  is  locally  observable  at  x^  under 
the  continuous  observatix?n  process  y(t)  for  all  t e Ct^,  t^D 
with  the  input  u. 

Proof:  From  Theorem  2.2  we  know  that  the  nonlinear  system 

is  locally  observable  under  the  discrete  observation  process 
Y at  Xq  with  input  u(t).  Therefore,  there  is  a one-to-one 
correspondence  between  Y (x^)  and  e 6 where  6 is  an  open 
neighborhood  of  x^.  Now  suppose  the  system  is  not  observable 
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under  the  continuous  observation  process  y(t)  = h(t,  x(t)); 

^ ^ Then  there  is  an  and  x ^ with  x^^, 

or  o o ' o 

2 19 

Xq  e 6 and  x x % such  that 
o o 

h(t,  g(t,  X u))  = h(t,  g(t,  x u) ) 

V O 

for  all  t e C^q/  this  contradicts  local  observ- 

ability under  the  discrete  observation  process  Y since  it 
implies 

Q.E.D. 

Similarly  we  can  make  the  following  statement  from 
Theorem  2.3. 

Theorem  2.5 

Let  y(t)  = h(t,  x(t))  be  the  continuous  observation  proc- 
ess associated  with  the  nonlinear  process 

x(t)  = f(t,  x(t),  u(t))  x(t  ) = X 

o o 

for  all  t e [ t^,  tj]. 

Let  {y(t^)},  i=l,  2,  t • »kbe  the  sequence  of  obser- 
vations at  any  k time  points  on  [t„,  t^l,  such  that  k > - 
and  y(ti)  has  a strong  derivative  at  Xq»  Then  the  non- 
linear system  is  locally  observable  at  x^  under  the  con- 
tinuous observation  process  h(t,  x(t));  t e [t^,  t^]  with 


Input  u if 


3y(ti) 


3y(t.) 


is  nonsingular. 

Application  of  the  foregoing  theorems  require  a mini- 
mum of  k = ^ observation  times  and,  in  general,  an  a priori 
knowledge  of  the  solution  so  that  an  expression  can  be 
obtained  for  y(tjj^)  = h(t^,  g(tj^,  x^,  u) ) . An  alternative 
test  for  local  observability,  which  avoids  these  problems, 
can  be  obtained  by  defining  a recursion  relationship.  Admis- 
sibility conditions  of  the  functions  f(*,*,*),  h {•,•),  and 
u(*)  as  defined  in  Section  2.1  will  have  to  be  made  somewhat 
more  restrictive.  For  the  following  theorems  f (•,*,•)  is 
continuous  and  has  continuous  partial  derivatives  with  respect 
to  all  of  its  arguments  up  to  and  including  order  n - 1 and 
continuous  mixed  partial  derivatives  of  the  same  order.  Also 
let  u(»)  be  continuous  and  have  continuous  derivatives  of 
order  n - 1.  Let  h(»,»)  be  continuous  and  have  continuous 
partial  derivatives  with  respect  to  both  of  its  arguments  up 
to  and  including  order  n and  continuous  mixed  partial  deriv- 
atives of  the  same  order.  Now  define  the  following 
recursion  relationship; 


y(t,  x)  = Fo(t,  x)  = h(t,  x(t)) 


^ = F,  (t,  X,  u)  = + f5_9,  f (t,  X,  u) 


d^y  3F,  fSF,')'*’  , 

^^2  = *.  «.  4)  = 7t"  + [sT-J  + [sitJ  ^ 


d^y 

dt” 


(2.1.14) 


An  nm  - vector  F(t,x)  is  formed  as  follows: 
F(t,x)  = [Fj^(t,x)Fj^^(t,x)  • • • ^n_:^(t,x)3 


(2.1.15) 


Theorem  2.6 

Let  Y be  the  discrete  observation  process  formed  by 
the  components  y(tj^),  i = 1,  2,  • • • Ic  with  tj^  e [t^,  t^^ 
which  is  associated  with  the  nonlinear  system 

x(t)  = f(t,  x(t),  u(t))  ~ ^o 

with  the  unique  solution  g(t,  x^,  u)  for  all  t e ]. 

Let  F(t,  x)  be  the  nm  - vector  defined  by  Eq  (2.1.15)  with 
components  obtained  from  the  recursion  relation  Eq  (2.1.14). 
If  for  one  of  the  observation  times  (say  t*  e {t^}) , F(t,  x) 
has  a strong  derivative  at  x(t*)  and  the  ran)c  of  the 
Jacobian  3F/8x  evaluated  at  t*  (denoted  by  F' (t*,  x(t*))  is 
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■s 

■if 


equal  to  n,  then  the  nonlinear  system  is  locally  observable 
under  the  discrete  observation  process  Y with  input  u. 


Proof:  From  the  definition  of  F,  we  know  that  F : x 

j^n  ^ r”™.  Also  from  Theorem  2.1,  if  F' (t*,  x(t*))  has  a 
rank  of  n,  then  an  inverse  function  exists  which  maps  points 
in  a neighborhood  of  F(t*,  x(t*))  to  a neighborhood  of 
x(t*);  i.e.,  F(t,  x)  is  locally  one-to-one  at  x(t*).  Now 
we  also  know  from  the  uniqueness  and  continuity  of  the 
solutions  x(t)  = g(t,  x^,  u)  that  there  is  a one-to-one 
mapping  from  a set  of  initial  conditions  in  an  open  neigh- 
borhood of  Xq*  which  yield  solutions  in  the  open  neighbor- 
hood 6*  of  x(t*).  Thus  for  observations  in  a neighborhood 
of  t*  we  have  established  the  local  one-to-one  mappings: 


relation  to  obtain  F(t*,  x(t*))  which  has  dimension  nm  and 


has  a rank  of  n.  Knowing  F(t*,  x(t*))  uniquely  determines 


x(t*)  and  x(t*)  uniquely  determines  x * 


— 

Note  that  this  theorem  relies  on  the  system  differential 
equation  having  a unique  solution  on  Ct^,  Appendix  B 

gives  a brief  summary  of  conditions  which  will  assure  a 
unique  solution. 

The  previous  theorem  can  be  easily  extended  to  a 
continuous  observation  process. 

Theorem  2.7 

Let  y(t)  = h(t,  x(t))  be  the  continuous  observation 
process  associated  with  the  nonlinear  system 

x(t)  = f(t,  x(t),  u(t))  ” ^o'  ^ ^ tj] 

with  a unique  solution  g(t,  x^,  u)  for  all  t e [t^,  t^]. 

If  there  exists  a t*  e [t  , t_]  such  that  F(t,  x)  has  a 

o r 

strong  derivative  at  x(t*)  and  the  rank  of  the  Jacobian 
3P/3x  evaluated  at  t*  is  equal  to  n,  then  the  nonlinear 
system  is  locally  observable  under  the  observation  process 
y(t),  for  all  t e L'^o'  ^f^  with  input  u. 

Proof:  Combine  Theorem  2.6  with  the  proof  of  Theorem  2.4. 

Consider  the  previous  ex^ple  where 


X^(t) 

z= 

X 

ft 

Xj(0)  = 

with 

2x2 (t) 

X2(0)  = x^2 

for  t e 

[0. 

• 

1 1 

The  observations  are 

= 

+ X2^(tj^)  i = 1,  2,  3. 
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Using  the  recursion  relation  Eq  (2.1.14),  we  get 


**0  *=  y(t,  x)  = + X2^(t) 

^1  = it  = + 2x2(t)X2(t) 

= 2x^2 (t)  + 4x2^(t) 

This  yields 


F(t,  x(t)) 


■ (t)  + x^^Ct) ■ 

2Xj^^(t)  + 4x2^(t) 


so  that 

F' (t,  x(t)) 


2xj^(t) 

4x^ (t) 


2X2 (t) 

8x^  (t) 


From  the  above,  we  see  that  the  F' (t,  x)  has  rank  two  for 
Xi(to)  ¥ 0 and  ^2 (t^)  ^ 0.  Thus  the  system  is  locally 
identifiable  for  x^  e {R^  ~ ^^ol  “ ® ^o2  ~ 8^}* 

2.3  Global  Observability  of  Nonlinear  Dynamical  Systems 

In  the  previous  section  we  discussed  local  one-to-one 
correspondence  and  observability.  Sufficient  conditions 
were  found  to  assure  a local  one-to-one  correspondence 
between  initial  conditions  in  a neighborhood  of  x^  and  the 
observation  process.  The  question  arises  as  to  what  con- 
ditions are  necessary  to  assure  that  a one-to-one  corre- 
spondence exists  for  any  x e X Cr”.  It  turns  out  that  a 


local  one-to-one  correspondence  at  each  point  x e x is  not 

o o 

sufficient.  A counter-example  to  such  a proposition  is 
obtained  from  the  previous  example.  That  is, 


Xj^(t)  = Xj^Ct) 


= 2x2(t) 


[0, 


(2.3.1) 


y(t^)  = 


i = 1,  2,  3 (2.3.2) 


W©  found  the  system  to  be  locally  observable  for 
2 

Xq  e {r  ~ {Xqj^  = 0 or  x^2  ~ 0}}.  Recall  that  the  solution 
has  the  form 


X (t)  = X . e^ 
1 ol 


' *02 


(2.3.3) 


It  is  clear  that  solutions  with  initial  conditions 


will  result  in  identical  observations;  i.e.,  Y(x  ^)  = Y (x  ^). 

o o 

This  exemplifies  the  need  for  additional  theory  to 
establish  conditions  for  global  observability. 


Definition;  The  mapping  F ; D Cr"  -v  r"  is  norm-coercive 
on  an  open  set  C D if,  for  any  Y > 0,  there  is  a closed, 
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Lil'l!! 


bounded  set  D c D such  that  | |F(x) | | > Y ^or  all 
Y o 

X e ~ D . 
o Y 

Note  that  if  = r”,  then  F is  norm-coercive  if  and 
only  if  lim  | |f(x) | | = “.  It  is  apparent  that  if 

lull  ^ ~ 

lim  I |f(x) II  = then  for  a y > 0 we  can  choose  a such 


that  for  X e R " I Ul  I will  be  sufficiently  large  so 
that  I |P(x) I I > Y«  On  the  other  hand,  if  F is  norm- 
coercive,  then  for  Y ■*■  “»  it  is  necessary  that 
lim  I |F(x)  II  = «>  to  ensure  that  | |f(x)  | | > Y for  all 

lull  ^ ” 

X e R*'  ~ D . 


Definition ; The  set  D is  path-connected  if,  for  any  two 
points  X,  y e D,  there  is  a continuous  mapping  p : [[0,l]]  D 

such  that  p(0)  = x and  p(l)  = y. 


Theorem  2 . 8 

Let  DCR*'  be  open  and  path-connected  and  assume  that 
F : DCR*^  -»•  r”  is  a local  homeomorphism  at  each  point  of  D. 
Then  F is  a global  homeomorphism  of  D onto  R^  if  and  only 
if  F is  norm-coercive  on  D. 


Proof;  See  Ref  39:137. 


Theorem  2.8  is  developed  based  on  the  continuation 
property  which  is  a relatively  standard  tool  in  the  theory 
of  ordinary  differential  equations. 
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Definition;  The  mapping  F ; DCR  -*■  r"  has  the  continuation 

property  for  a given  continuous  function  q : CR^  -►  r” 

if  the  existence  of  a continuous  function  p ; Co,a)  -►  D, 

a e (0,1]  such  that  F(p(t))  = q(t)  for  all  t e []0,a)  implies 

that  lim  _p(t)  = p(a)  exists  with  p(a)  e D,  and  F (p (a) ) =q (a) . 
t-»-a 

A discussion  of  this  property  is  contained  in  Ref  39:133. 

In  all  of  the  theorems  where  norm-coerciveness  is  stated, 
one  can  replace  "norm-coercive"  by  "has  the  continuation 
property"  and  the  theorem  will  remain  valid.  Norm-coer- 
civeness  is  stated  in  the  theorems  since,  in  most  cases,  the 
norm-coercive  condition  is  probably  easier  to  verify  than 
the  continuation  property.  It  can  be  shown  that  the  con- 
tinuation property  holds  for  all  continuously  differentiable 
functions  (Ref  39:140).  Recall  from  above  that  "continu- 
ously differentiable"  means  a continuous  Gateaux-derivative 
(and  thus  a continuous  Frechet  derivative) . Combining  this 
with  Theorem  2.8  results  in  the  following  useful  theorem. 

Theorem  2.9 

Let  D be  open  and  path-connected  and  assume  that 
F ; DCR*^  ->  r”  is  continuously  differentiable  on  D.  If  the 
Jacobian  matrix  F' (x)  is  nonsingular  for  all  x e D,  then  F 
is  one-to-one  on  D. 

With  the  above  definitions  and  theorems  as  background, 
we  can  now  make  statements  similar  to  Theorems  2.8  and  2.9 
for  a mapping  g : X Cr”  Y CR^  where  k > n. 
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Theorem  2.10 

Suppose  g : X C r”  ->  where  k > n.  Let  X C X 

- pc 

where  X C r"  is  open  and  path-connected.  Assume  that  for 
pc 

all  X e Xp^  g has  a strong  derivative  at  x and  we  can  find 
a common  n x n submatrix  of  dg/dx  which  is  nonsingular. 

(By  common  submatrix,  we  mean  it  is  always  composed  of  the 
same  rows  of  dg/dx)  . Define  a mapping  g^^  : Xp^  C r”  -»■  r” 
where  g^^  (x)  is  composed  of  the  n components  of  g (x)  which 
correspond  to  the  nonsingular  n x n submatrix  of  dg/dx. 

Then  g is  a one-to-one  mapping  of  X onto  R^  if  g^^  is  norm- 
coercive  on  Xp^. 


Proof:  If  k = n,  then  we  let  X = Xp^  and  g = g^^  so  that 

this  becomes  a statement  of  Theorem  2.8. 

For  k > n,  suppose  g^  is  norm-coercive  on  Xp^.  Then 

from  Theorem  2.6  we  know  g^^  is  a homeomorphism  of  Xp^  onto 

r’’.  Now  suppose  g is  not  a one-to-one  mapping  of  X onto 
k 

R . Then  there  exist  Xj^,  X2  e X with  x^^  / X2  such  that 
g(x^)  = g(x2).  But  this  contradicts  g^^  being  a homeomorphism 
since  from  the  definition  of  g^^,  this  would  require 
gi(xi)  = gi(x2) . 

Q • E • D • 

Note  that  if  X = the  norm-coerciveness  condition 

can  be  replaced  by  lim  1 |g2(x)  | I = «>  in  Theorem  2.10. 

11x11  00 

Corollary:  Suppose  g : X C r”  ->•  R^  where  k ^ n and  X is 

open  and  path-connected.  If  there  exists  a g^^  formed  from 
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n components  of  g such  that  is  continuously  differentiable 
on  X and  g^^*  (x)  is  nonsingular  for  all  x e X,  then  g is  a 


one-to-one  mapping  on  X. 

Proof:  Follows  directly  from  Theorems  2.9  and  2.10. 

Theorem  2.10  and  its  corollary  allow  us  to  state  suffi- 
cient conditions  for  a nonlinear  process  to  be  globally 
observable  under  the  observation  process  Y for  all  Xn  £ X^ 
with  input  u.  We  can  make  such  a statement  when  there  is  a 
one-to-one  correspondence  between  any  e Xq  C r”  and  the 
observation  vector  Y for  an  input  u.  Refer  to  Section  2.1 
for  the  definitions  of  the  nonlinear  process  and  admissi- 


correspondence  between  Y and  any  initial  state  vector 


^0^*0^  *pc  norm-coercive  on  Xpc  or  (2) 

Yj^  is  continuously  differentiable  on  Xp^. 


Proof:  This  follows  directly  from  Theorem  2.10  and  its 

corollary. 


Note  that  if  in  Theorem  2.11,  then  the  norm- 

pc  ' 

coercive  condition  is  replaced  by 

lim  I |Yi(Xq)  I I = «» 

I IxqI  I “ 

Theorem  2.11  gives  sufficient  conditions  for  a nonlinear 
process  to  be  globally  observable.  For  some  problems,  a 
more  practical  set  of  conditions  can  be  obtained  from  the 
properties  of  monotone  functions. 

Definition:  A mapping  F : D C r”  -»•  r”  is  monotone  on 

Do  C D if 

[f(x)  - FCy)]"^  [x  - y]  ^ 0 

for  all  X,  y e Do»  F is  strictly  monotone  on  Do  if  strict 
inequality  holds  whenever  x y,  and  uniformly  monotone  if 
there  is  a Y > 0 so  that 

[f(x)  - F(y)]*^  [x  - y]  > Y(x  - y)^  (x  - y) 

for  all  X,  y e D^. 
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Theorem  2.12 

Let  F : D C r”  R”  be  strictly  monotone  on  Dq  C D. 

Then  P is  a one-to-one  mapping  on  D^. 

Proof:  From  the  definition  of  strict  monotonicity 

[p(x)  - F(y)]'*’[_x  - y3  > 0 

for  all  X,  y e Dq  X y. 

If  F were  not  one-to-one,  then  there  exists  an  x,  y e Dq 
such  that  F(x)  = F(y)  with  x y.  This  contradicts  the 
defining  inequality. 

Q.E.D. 

Theorem  2.13 

n ]c 

Suppose  g : X C R R where  k ^ n.  If  there  exists 
a g-^  formed  from  n components  of  g such  that  g^^  is  strictly 
monotone  on  X,  then  g is  one-to-one  on  X. 


components  of  Y such  that  is  strictly  monotone  on 

Xq  C R , then  there  is  a one-to-one  correspondence  between 

y and  any  initial  state  vector  x e X C R^. 

o o 

Proof:  Follows  directly  from  Theorem  2.13. 

Another  set  of  conditions,  which  may  be  of  use  in  some 
cases,  can  be  derived.  First  a couple  of  introductory 
theorems  will  be  established. 

Definition:  A subset  D C r”  is  convex  if,  whenever  it 

o 

contains  x and  y,  it  also  contains  Xx  + (1  - X)y  for  all 
X,  0 £ X £ 1.  Note  that  the  set  {z  : z = Xx  + (1  - X)y 
for  0 ^ X £ 1}  is  called  the  line  segment  joining  x and  y. 

Theorem  2 . 15 

Let  F : D C R*^  -»•  r”  be  continuously  differentiable  on 

an  open  convex  set  C D.  If  dF/dx  is  positive  definite 
for  all  X e Dq,  then  F is  strictly  monotone  on  D^. 

Proof:  See  Ref  39:142. 

Corollary:  IfF:DCR”  -^R^'is  continuously  differen- 

tiable on  an  open  convex  set  C D and  dF/dx  is  positive 

definite  for  all  x e D^,  then  F is  one-to-one  on  D^. 

o o 

Proof:  From  Theorem  2.15,  F is  strictly  monotone  on  D^. 

From  Theorem  2.12,  this  is  sufficient  to  assure  F is  one- 
to-one  . 
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Following  the  same  approach  that  was  used  In  arriving 
at  Theorem  2.14  from  Theorem  2.12,  we  obtain 

Theorem  2.16 

Let  Y be  the  mk  observation  vector  associated  with  the 
nonlinear  process 

x(t)  = f(t,  x(t),  u(t))  x(t^)  = Xq  e Xq 

so  that  Y : XqC  r”  -*■  Define  Y^  : C r"->-r”  by  taking 

any  n components  of  Y (x^) . Then  there  is  a one-to-one 
correspondence  between  Y and  any  initial  state  vector 
Xq  e Dq  C Xq  (i.e.,  the  nonlinear  system  is  globally  observ- 
able on  Dq  under  the  observation  process  Y for  input  u)  if 
Yl  is  continuously  differentiable  on  an  open,  convex  set 
Dq  C Xq  and  3Y3^/3xq  is  positive  definite  for  all  x^  g D^. 

Theorems  2.11,  2.14,  and  2.16  can  be  used  to  make 
statements  of  sufficient  conditions  for  global  observability 
of  nonlinear  systems  under  continuous  observations. 

Theorem  2.17 

y(t)  = h(t,  x(t))  be  the  continuous  observation 
process  associated  with  the  nonlinear  system 

x(t)  = f(t,  x(t),  u(t))  “ ^o  ^o 

for  all  t e [t^,  . 

Let  Y be  the  discrete  observation  process  vector  formed  by 
the  components  y(tj^),  i = 1,  2,  • • • k with  k > ^ and  t^ 

33 


any  k time  points  on  C^q»  Then  the  nonlinear  process 

is  globally  observable  on  C r”  under  the  continuous 

observation  y(t)  for  all  t e ft  , t_]  with  input  u,  if  one 

o f 

of  the  following  conditions  is  met; 

(1)  Xq  is  open  and  path-connected  and  there  exists  a 
formed  from  n components  of  Y such  that  has  a strong 

derivative  at  and  3Yj^/3x  is  nonsingular  for  all  x^  e X^ 
and  Yj^  is  norm-coercive  on  Xq. 

(2)  Xq  is  open  and  path-connected  and  there  exists  a 
Yj^  formed  from  n components  of  Y such  that  Yj^  is  continu- 
ously differentiable  on  X^  and  ^Yj^/3x^  is  nonsingular  for  all 

e X^. 

o o 


(3)  There  exists  a Yj^  formed  from  n components  of 

Y such  that  Y,  is  a strictly  monotone  function  of  x for  all 
-1-  o 


^o  ^ ^o* 


(4)  X^  is  an  open,  convex  set  and  there  exists  a Yj^ 
formed  from  n components  of  Y such  that  Yj^  is  continuously 
differentiable  and  9^i/9Xq  is  positive  definite  for  all 


Proof;  Follows  directly  from  Theorems  2.11,  2.14,  and 
2.16  using  the  same  proof  as  used  for  Theorem  2.4. 


Using  the  recursion  relation  of  Eq  (2.1.14),  we  can 
formulate  another  set  of  sufficient  conditions  for  global 
observability. 


Theorem  2.18 


r 


Let  Y be  the  discrete  observation  process  formed  by 
the  components  y(t^),  i = 1,  2,  • • • k which  is  associated 
with  the  nonlinear  system 


x(t)  = f (t,  x(t) , u(t) ) 


x(t  ) = X e X 
o o o 


with  a unique  solution  g(t,  x , u)  for  all  t e(^t  , t ]. 

o of 


Let  F(t,  x)  be  the  nm-vector  defined  by  Eq  (2.1.15).  Let 
t*  denote  any  one  of  the  observation  times  {tj^}  and  let  the 
solution  g(t,  u)  map  the  initial  conditions  x^  e onto 

X*  at  time  t*.  Then  the  nonlinear  process  is  globally 
observable  on  X^  C R*'  under  the  observation  process  Y with 
input  u if  one  of  the  following  conditions  is  satisfied: 

(1)  X*  is  open  and  path-connected  and  there  exists  an 
Fj^(t,  x)  formed  from  n components  of  F such  that 

Fj' (t*,  x(t*))  is  strong  and  is  nonsingular  for  all 
x(t*)  e X*  and  Fj^(t*,  x(t*))  is  norm-coercive  on  X*. 

(2)  X*  is  open  and  path-connected  and  there  exists  an 
F^(t,  x)  formed  from  n components  of  F(t,  x)  which  at  t*  is 
continuously  differentiable  on  X*  and  Fj^' (t*,  x(t*))  is 
nonsingular  for  all  x(t*)  e X*. 

(3)  There  exists  an  F^^  (t,  x)  formed  from  n components 
of  F(t,  x)  such  that  F(t*,  x(t*))  is  a strictly  monotone 
function  of  x(t*)  for  all  x(t*)  e X*. 

(4)  X*  is  an  open,  convex  set  and  there  exists  an 
Fj^Ct,  x)  formed  from  n components  of  F(t,  x)  such  that 
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x(t*))  is  continuously  differentiable  and 
(t*,  x(t*))  is  positive  definite  for  all  x(t*)  e X*. 

Proof:  Combine  Theorems  2.11,  2.14,  and  2.16  with  the 

proof  of  Theorem  2.6. 

Theorem  2.18  can  be  extended  to  a continuous  obser- 
vation process  very  easily. 

Theorem  2 . 19 

Let  y(t)  = h(t,  x(t))  for  all  t e be  the 

continuous  observation  process  associated  with  the  nonlinear 
system 

x(t)  = f(t,  x(t),  u(t))  SfCt^)  = Xq  ^ ^o 

for  all  t e t£3« 

Let  t*  e ^f^  9(t*»  u)  map  the  initial  con- 

ditions Xq  e Xjj  onto  X*  at  t*.  Then  the  nonlinear  process 
is  globally  observable  on  Xq  C under  the  continuous 
observation  process  y(t)  with  input  u if  any  one  of  the 
conditions  of  Theorem  2.18  is  satisfied. 

Proof:  Combine  Theorem  2.18  with  the  proof  of  Theorem  2.4. 

Consider  the  example  presented  at  the  beginning  of  this 
section. 


y(t^)  = X it^)  + X (tj^)  i = 1,  2,  3 

so  that  the  observation  process  is 


2 2 . 

Y = 

Xi  (t^)  + x^  (t^) 
2 2 , . 

Xi  (t3>  + x^  (t3)^ 

2 2ti  2 4ti 

X,  e ■'•  + x_  e 
ol  o2 

= 

Xo2^  ^2^2  + 

g4t2 

2 2t,  2 

4t, 

X - e •^  + X ^ 
ol  o2 

5 

We  showed  that  the  nonlinear  system  was  locally  observable 

2 

under  the  observation  process  on  = {R  ~ ^^ol  ~ ® 

Xq2  = 0)};  however,  the  system  is  not  globally  observable 

on  X^.  Suppose  we  test  for  global  observability  on  the  set 
2 

Xq'  C R where  X^'  = {x^  d x^^j^  > 0,  x^2  > 0^*  That  is 

2 

X^'  is  the  following  subset  of  R : 


system  is  globally  observable  if  Y-i  is  norm-coercive  on  X ' . 

o 

If  is  norm- coercive,  then  for  any  Y > 0,  there  must  be  a 
closed  bounded  set  C X^'  such  that  1 | | > Y for 

Xq  e X ' - D^.  From  Y^  we  see  that  | |Yj^|  | > e^^l  for 

Xq  e Xq,  so  that 


Therefore  we  want  to  find  an  x.,(y)  such  that 


From  this  we  see  that  given  a y > 0,  we  choose 


^ = t^o'  ~ t-oi  > Y;%  ^ ^02  > o» 


This  shows  is  norm-coercive  on  X^*  and  therefore  the 
system  is  globally  observable  on  Xq*.  This  illustrates  that 
for  nonlinear  processes,  there  may  be  observability  on  only 
parts  of  r”.  Also,  although  not  shown  by  the  example,  it  is 
clear  that  observability  of  a nonlinear  system  depends  on 
the  input  function  u. 


2,4  Observability  of  Linear  Time  Invariant  Systems 
Consider  systems  which  are  described  by 


x(t)  = A x(t)  + B u(t) 


and  an  observation  process 


t E [t^,  t^],  x(t^)  = X^ 


(2.4.1) 


y(t^)  = C x(t^)  + D u(t^)  i = 1,  2, 


(2.4.2) 


Suppose  we  apply  Theorem.  2.6  to  this  problem  to  test  for 
observability.  From,  the  recursion  relation  Eq  (2.1.14): 


= y(t)  = C x(t)  + D u(t) 

= C x(t)  + D u(t) 

= C A x(t)  + C B u(t)  + D u(t) 

= C A x(t)  + C B u(t)  + D*u(t) 

* C A^  x(t)  + C AB  u(t)  + C B u(t)  + D*u(t) 


= C x(t)  + C A^B  u(t)  + C AB  u(t)  + C B li(t)  + D u(t) 


F = C a""^  X(t)  + C a"“^B  u(t)  + C A^-^B  u(t) 
n— 1 


3”"^2 


+ C B 


3tn-2  + D g^n-1 


Also  recall  that 


= [ 


X,  T „ T 
^1 


n-1  J 


and  therefore 


3F 


3F„  I 3Fj^  I 


I 3P 


3x  I 


I • 


3x 


r T I T T I T~  T ! ' , T1 

II  = [c  ; A C 1 A C , • . • ; A ‘^C  J 


(2.4.3) 


Thus  from  Theorem  2.6,  the  system  is  locally  observable  if 
3F/3x  has  rank  n. 

To  test  for  global  observability  we  can  use  Theorem 
2.18.  In  addition  to  the  condition  3F/3x  being  of  rank  n, 
we  must  show  that  any  set  of  n components  of  F is  norm- 
coercive  on  the  solution  set  X.  For  this  case,  the  obser- 


vation process  is  a mapping  Y(Xq)  : x^  e r”  -»•  Recall 


n 


from  page  26  that  for  this  case  F is  norm-coercive  on  R if 


lim  I |F(t,  x(t) ) I I = “> 
I |x  (t)  I I 
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This  IS  seen  to  be  true  by  looking  at  the  components 

previous  page.  From  the  above  we  see 
that  linear  time  invariant  systems  are  globally  observable 
on  R if 


m I T T * T * 

C C I AC  I C I 
I I I 


is  of  rank  n.  This  result  agrees  with  that  obtained  from 
linear  system  theory.  Note  that  observability  does  not 
depend  on  the  input  function  u. 

Clearly  these  same  results  apply  for  a continuous 
observation  process. 


2.5  Summ.ary 

In  this  chapter  sufficient  conditions  have  been  estab- 
lished for  determining  local  and  global  observability  of 
nonlinear  systems. 

Local  Observability.  Theorems  2.2  and  2.3  state 
sufficient  conditions  for  local  observability  of  nonlinear 
system.s  under  discrete  time  measurements.  Theorem  2.6 
provides  alternative  conditions  based  on  a recursion  rela- 
tionship which  can  avoid  the  requirement  of  an  a priori 
knowledge  of  the  solution,  but  requires  the  system  differ- 
ential equation  have  a unique  solution  on  Corre- 

sponding sufficient  conditions  for  establishing  local 
observability  under  continuous  time  measurements  are  stated 
in  Theorems  2.4,  2.5,  and  2.7. 


Global  Observability.  Theorem  2.11  extends  the  results 


of  Theorem  2.3  for  discrete  time  measurements  to  sufficient 
conditions  for  global  observability  based  on  a norm-coercive 
or  a continuously  differentiable  condition.  Theorems  2.14 
and  2.16  state  sufficient  conditions  for  global  observ- 
ability under  discrete  time  measurements  using  the  properties 
of  strictly  monotone  and  convex  functions  respectively. 
Theorem  2.18  extends  Theorem  2.6  to  sufficient  conditions 

i 

for  global  observability  under  discrete  measurements  based 
on  a recursion  relation  by  using  the  results  of  Theorems 
2.11,  2.14,  and  2.16.  Finally  Theorems  2.17  and  2.19  extend 
the  above  results  to  sufficient  conditions  for  global  observ- 
ability under  continuous  time  measurements. 


Chapter  III 


IDENTIFIABILITY  OF  NONLINEAR  DYNAMICAL  SYSTEMS 


The  previous  chapter  developed  new  results  concerning 
observability  of  nonlinear  systems  under  discrete  and  con- 
tinuous observations.  We  now  want  to  apply  these  results 
to  the  problem  of  identifying  parameters  of  nonlinear 
dynamical  systems. 

3.1  Nonlinear  State  Equations  and  Definition  of 
Identif iability 

Consider  the  nonlinear  systems  described  by 

> 

* x(t)  = f(t,  x(t),  u(t),  0)  x(t^)  = Xq  (3.1.1) 

and  the  discrete  observation  process 

y(tj^)  = h(ti»  x(ti),  4>)  i = 1,  2,  • • • , k 

(3.1.2) 

where  x(t)  is  an  n vector,  y(tjj^)  is  an  m vector,  u is  an  r 
vector  describing  the  input  function,  is  an  s vector  of 
constant  parameters,  and  t^  e [t^,  t^]  for  all  i.  Assume 
£(•»•»•»•)  has  continuous  partial  derivatives  with  respect 
to  t,  X,  and  <|>.  u is  a continuous  function  on  Ct^,  t^D 
and  h (•,»,•)  has  continuous  partial  derivatives  with  respect 
to  its  three  arguments. 

A solution  of  Eg  (3.1.1)  is  denoted  by  g(t,  x^,  u,  <J>) . 

In  the  identif iability  problem,  we  are  interested  in  finding 
sufficient  conditions  to  guarantee  a one-to-one  correspondence 
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between  the  parameter  vector  (fi  and  the  observed  sequence 
{y(t^)}.  That  is,  we  want  to  assure  that  the  value  of  (p 
is  uniquely  determined  from  the  knowledge  of  u(t)  and  the 
output  sequence  {y(t^)}  for  all  t e [t^,  t^] . As  in 

Chapter  II,  the  observation  process  is  viewed  as  a single 
mk  vector 


Y = Cy’’(t,)y'^(t2)  . . . 


(3.1.3) 


The  parameter  vector  <|)  of  the  nonlinear  system  is 
identifiable  under  the  observation  process  Y with  input  u 
if  there  is  a one-to-one  correspondence  between  the  param- 
eter vector  ({)  and  the  observation  vector  Y(4>).  The  param- 
eter vector  is  locally  identifiable  at  <l>  when  the  one-to- 
one  correspondence  holds  in  a neighborhood  of  <|) . The 
parameter  vector  is  globally  identifiable  on  4*^  under  the 
observation  process  Y when  there  is  a one-to-one  corre- 
spondence between  (p  and  Y for  all  (Ji  e C R®  with  the 
input  u.  By  defining  an  augmented  t'tate  vector 


Xa(t) 


x(t)* 


(3.1.4) 


it  becomes  clear  that  identification  of  <|)  is  a special  case 
of  the  observability  problem.  That  is 


Xa(t) 


f (t,  x(t)  , u(t)  , (j))“] 


(3.1.5) 
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* 


with 


*a<*o>  = 


x(t^) 


(3.1.6) 


Thus  identifying  <l>  becomes  part  of  the  observability  problem 
for  X (t) . In  this  chapter,  we  will  assume  the  initial  con- 

Si 

ditions  x(t^)  are  known  and  use  the  results  of  Chapter  II  to 
state  conditions  for  identif lability . 


3.2  Local  Identifiab i lity  of  Non linear  Dynamical  Systems 


In  this  section,  we  will  follow  the  same  approach  as 
Section  2.2  in  deriving  sufficient  conditions  for  the  param- 
eter vector  <}»  of  a nonlinear  system  to  be  locally  identi- 
fiable. 


Theorem  3 . 1 

Let  Y be  the  mk  observation  vector  associated  with  the 
nonlinear  process 


x(t)  = f(t,  x(t),  u(t),  (|)) 


with  the  dimension  of  Y equal  to  or  greater  than  the  param- 
eter vector  i.e.,  mk  ^ s and  let  Y have  a strong  deriv- 
ative at  Xq.  Then  the  parameter  vector  is  locally  observ- 
able under  the  observation  process  Y with  input  u and  ini- 
tial conditions  x^  if  3Y/3<|i_-hae-  rank  s. 


Proof;  This  follows  directly  from  Theorem  2.1, 


Theorem  3.2 


Let  Y have  a strong  derivative  at  x^.  Then  the  param- 
eter (|>  is  locally  identifiable  under  the  observation  proc- 
ess Y with  input  u and  initial  conditions  x if 


is  nonsingular. 

Proof;  This  is  a special  case  of  Theorem  2.3. 


These  theorems  can  be  easily  extended  to  continuous 
observations  as  follows; 


Theorem  3.3 

Let  y(t)  = h(t,  x(t),  4>)  for  all  t e [ t^,  t^3  he  the 
continuous  observation  process  associated  with  the  nonlinear 
process 

x(t)  = f(t,  x(t),  u(t),  (}>)  x(t  ) = X 

o o 

Let  Y be  the  discrete  observation  vector  formed  from 

i = 1,  2,  • • • , k,  the  sequence  of  observations  at 
any  k time  points  on  [t^,  t^]  such  that  k > and  let  Y have 
a strong  derivative  at  x^.  Then  the  parameter  (J>  is  locally 
identifiable  under  the  continuous  observation  process  y(») 
for  all  t E [t^/  tjl  with  input  u and  initial  conditions 


I 


i 


if  either  of  the  following-  conditions  is  met; 
(1)  3y/94>  has  rank  s. 


k 

(2)  I 

i=l 


■3y(tj.)' 

T 

■ay(ti)' 

34> 

34) 

is  nonsingular. 


Proof;  This  follows  from  Theorems  2.4  and  2.5. 


The  recursion  relation  of  Eq  (2.1.14)  can  be  used  to  derive 
additional  sufficient  conditions  for  local  identifiability . 
The  recursion  can  be  used  to  derive  the  functions  F , F,  , 

F2  • * • » ^n+s-1*  vector  is  formed; 

f(t,  4))  = [F^^(t,  <l))F3^'’’(t,  4,)  . . . Fn+s_i'^(t,  4>)]  . Fol- 

lowing the  same  approach  as  Theorem  2.6,  we  can  state  the 
follov/ing  theorem. 

Theorem.  3 . 4 

Let  Y be  the  discrete  observation  process  formed  by 
the  components  y(t^),  i = 1,  2,  • • • k with  t-^  e [t^,  t^] 
which  is  associated  with  the  nonlinear  system. 

x(t)  = f(t,  x(t),  u(t),  4»)  5c(t^)  = x^ 

with  a unique  solution  g(t,  x , u,  4>)  for  all  t e ft  , t J . 

o of 

Let  F(t,  4>^  be  the  (n  + s)m  vector  formed  from,  the  recursion 
relation  Eq  (2.1.14) . If  for  one  of  the  observation  times 
(say  t*  e {tj^}),  F(t*,  4>)  has  a strong  derivative  at  x^  and 

the  rank  of  the  Jacobian  3F/34>  evaluated  at  t*  is  equal  to 
s,  then  the  nonlinear  system  is  locally  identifiable  under 
the  discrete  observation  process  Y with  input  u and  initial 
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conditions  x^. 

Proof:  Follows  from  Theorem  2.6. 

A similar  statement  can  be  made  for  continuous  observations. 
Theorem  3.5 

Let  F(t,  (())  have  a strong  derivative  at  x^.  If  for  a 
continuous  observation  process,  the  rank  of  the  Jacobian 
9F/3({i  evaluated  at  t*  e CiQ/  t^]  is  equal  to  s,  then  the 
nonlinear  system  is  locally  identifiable  under  the  contin- 
uous observation  process  y(*)  for  all  t e C^o'  ^f^  with 
input  u and  initial  condition  x^. 

3.3  Global  Identif iability  of  Nonlinear  Dynamical  Systems 

As  with  the  observability  problem,  local  identifi- 

ability  at  each  point  ()>  e is  not  sufficient  to  assure  that 

a one-to-one  correspondence  exists  between  (j)  and  Y for  any 
s 

^ e $ C R . We  can  draw  upon  the  theorems  in  Chapter  II 
concerning  global  observability  to  state  some  theorems  on 
global  identif iability. 

Theorem  3 . 6 

Let  Y be  the  mk  observation  vector  associated  with  the 
nonlinear  system 

x(t)  = f(t,  x(t),  u(t),  (|>)  xCt^)  = Xq 

The  parameter  vector  4>  is  globally  identifiable  on  4>  under 
the  observation  process  Y with  input  u and  initial 
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conditions  if  one  of  the  following  conditions  is 
satisfied: 

(1)  $ is  open  and  path-connected  and  there  exists  a 
formed  from  s components  of  Y such  that  Yj^  has  a strong 

derivative  at  x^  and  8Yj^/3(|)  is  nonsingular  for  all  <|)  e $ 
and  Y^  is  norm-coercive  on  4>. 

(2)  $ is  open  and  path-connected  and  there  exists  a 
Yj  formed  from  s components  of  Y such  that  Y^  is  contin- 
uously differentiable  and  3Y2^/3(|)  is  nonsingular  for  all 

<t>  e 

(3)  There  exists  a Yj^  formed  from  s components  of  Y 
such  that  Yi  is  a strictly  monotone  function  of  (()  for  all 
4)  e 4>. 

(4)  4»  is  an  open,  convex  set  and  there  exists  a Yj^ 
formed  from  s components  of  Y such  that  Y]l  continuously 
differentiable  and  3Yj^/3(|)  is  positive  definite  for  all 

4>  c 4. 

Proof:  Follows  from  Theorems  2.11,  2.14,  and  2.16. 

Following  the  Scune  approach  as  for  the  observation 
problem,  we  can  extend  the  above  theorem  to  the  continuous 
observation  process  y(t)  = h(t,  x(t))  for  all  t e [t^,  tj] . 
Theorem  3,6  applies  to  this  case  by  letting  Y be  formed 
from  any  k observations  times  of  the  continuous  observation 
on  [tj,,  tj]  with  k > |. 
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The  recursion  relation  of  Eq  (2.1.14)  can  be  applied 
to  the  global  identification  process.  As  with  the  local 
identification,  the  recursion  relation  can  be  used  to 
derive  the  functions  F2,  • • • ^n+s-1*  Then  an 

(n  + s)m  vector  is  formed  as  follows; 

F(t,  <|))  = (|))Fi'^(t,  <!>)•••  (|))]'^ 

(3.3.1) 


Theorem  3 . 7 

Let  y be  the  discrete  observation  process  formed  by 
the  components  y(t^),  i = 1,  2,  • • , k with  t^^  e 

which  is  associated  with  the  nonlinear  system 

x(t)  = f(t,  x(t),  u(t),  4))  “ *0 


with  a unique  solution  g(t,  x , u,  A)  for  all  t e [t  , t_l. 

o o 

Let  F(t,  <|))  be  the  (n  + s)m  vector  formed  from  the  recursion 

relation  as  indicated  in  Eq  (3.3.1).  Let  t*  denote  any  one 

of  the  observation  times  {tj^}  so  that  F(t*,  4>)  is  the 
• Z * ^ ..s  „(n+s)m 

mapping  F:$CR  R « Then  the  nonlinear  process 

is  globally  identifiable  on  4>  under  the  observation  process 
Y with  input  u and  initial  conditions  x^  if  one  of  the  fol- 
lowing conditions  is  satisfied; 

(1)  4>  is  open  and  path-connected  and  there  exists  an 

F^(t,  (^)  formed  from  s components  of  F(t,  (j))  such  that 
3F/3(|i  evaluated  at  t*  is  strong  and  is  nonsingular  for  all 
♦ e ♦ and  Fj^(t*,  <^)  is  norm-coercive  on  ♦. 


(2)  t is  open  and  path-connected  and  there  exists  an 
Pj^Ct,  ♦)  formed  from  s components  of  F(t,  4>)  which  at  t*  is 
continuously  differentieible  and  3Fj^/34)  is  nonsingular  for 
all  e ». 

(3)  There  exists  an  F^^Ct,  <t>)  formed  from  s components 

of  F(t,  4»)  such  that  F(t*,  <|))  is  a strictly  monotone 
function  for  all  e 

(4)  4 is  an  open,  convex  set  and  there  exists  an 
Fj^(t,  (>)  formed  from  s components  of  F(t,  (}>)  such  that  Fj^ 
is  continuously  differentiable  and  3Fq^/3(|)  evaluated  at  t* 
is  positive  definite  for  all  (pet. 

Proof:  Follows  from  Theorem  2.18. 

For  continuous  observations,  the  following  theorem  can 
be  stated. 

Theorem  3.8 

Let  y(t)  = h(t,  x(t),  4>)  for  all  t e t^]  be  the 

continuous  observation  process  associated  with  the  non- 
linear system 

X(t)  = f(t,  X(t),  U(t),  (J))  = Xq 

with  a unique  solution  g(t,  x , u,  <|>)  for  all  t e ft  , t,]. 

o of 

Let  t*  c Ci-Q/  tfl  so  that  F(t*,  4>)  is  the  mapping 
P ; ♦ C R®  Then  the  nonlinear  process  is  glob- 

ally identifiable  on  t under  the  continuous  observation 


process  y(t)  for  a2J.  t c [t^,  t^]  with  input  u and  initial 
conditions  Xq  if  any  one  of  the  conditions  of  Theorem  3.9  is 
satisfied. 

3.4  Identifiability  of  Linear  Time  Invariant  Systems 

Consider  systems  which  are  described  for  all  t e [t^,  t^] 

x(t)  = A((|))x(t)  + B((j))u(t)  x(t  ) = x^  (3.4.1) 

o o 


with  a discrete  observation  process 
y(ti)  = C(<{))x(ti)  + D(<{))u(ti)  i = 1,  2,  • • 'k  (3.4.2) 

Suppose  we  use  Theorem  3.4  to  test  for  local  identifiability 
of  (|).  From  the  recursion  relation 


Fq  = C(4))x(t)  + D(<|))u(t) 

F^  = C(4))x(t)  + D((|))u(t) 

= C((|))A((|))x(t)  + C((^))B((^)u(t)  + D(<J))u(t) 

F2  = C((J>)A((t))x(t)  + C((|))B((|))u(t)  + D((|))u(t) 

= C((l))A2((J))x(t)  + C(4))A(4))B((|))u(t)  + C((J))B((t))u(t) 

+ D(<|))u(t) 

F3  = C(<|))A^  (<J))x(t)  + C(<j))A^((())B((|))u(t)  + C((J))A(4i)B((())u(t) 

+ C(<J))B(()))u(t)  + D(()))*u*(t) 


p 

'^n+s-l 


= C((J))A"'^®“^(().)x(t)  + C((|))A"'^®“^(<|))B(()>)u(t) 


n+s-3 


+ C(4>)a"^®  •^B(4))u(t)-  • • + C((|))A((|))B(<)))£ H 


~n+s-2„  ~n+S“l 

+ C(4>)B((J.)5 — — + D(().)1  " 


3 t 


n+s-3 


3 ^n+s-2 
52 


3 tn+s-1 


From  this  we  obtain 


.i 


p 


C(<j))x(t)  + D((l))u(t) 

C((|>)A((|))x(t)  + C(<J))B(((.)u(t)  + D(<}))u(t) 
C(«|»)A^(()))x(t)  + C(<|))A((|))B((I))u(t)  + C((J))B(<())u 
+ D(4>)u(t) 


,n+s-l 


((|))x(t)  + C(4))A 


n+s-2 


(4))B((t))u(t)  • 


+ D((|)) 


3nts-l^(t) 

9tn+s-l 


C((j)) 

C(4))A((|)) 

C(<|.)a2((|)) 

• n+s-1 
C(<^)A  (4>) 


D((|)) 

C(4))B((|)) 

C(c{>)A((|))B((|)) 

• n+s-2 

C((J))A  ((t)B((l)) 


0 

D(<J>) 

C(<|))B((|)) 

• n+s-3 

C(<|))A  ((())B(4)) 


Note  that  the  nonzero  elements  of  through  F^^^g  are 
contained  in  F2. 

From  above  we  have 

.n+s-*!  ... 

_ _ - - 8 u(t) 

F = F^  x(t)  + Fj  x(t)  + • • • + F^^g  g^n+s-I 


Thus 


the  derivative  of  F with  respect  to  the  vector  (!>  is 


1 3Fi 
1 ^^2 

1 

x(t)  .1  • • • 

1 

1 3Fi 
1 ^ 

x(t)| 

\ 

3F,  1 

t 

, ,1 

u(t)  1 
1 

' • 

3<P2  I 

• • 1 

1 

551  ““’J  ■ 

S^n+s 

1 ’'■n.s 

3<I>1 

g^n+S-1  1 

1 at"'*'®" 

^u(t) 

-1 


From  Eq  (3.4.3),  it  is  seen  that  3F/3<|)  will  have  rank  s if 
either  of  the  following  conditions  is  satisfied: 

3F, 

(1)  has  rank  s 


(2)  has  rank  s 


where 


^1' 


C(<|.) 

C((J))A(())) 
C((|))a2  {())) 


C((j))A"'*'®"^(<l))J 


’ 


■d((|)) 

C(((.)B((j)) 

C((())A(<j))B((j)) 

.C(4.)A”‘^®“^(<|))B((t)). 


and 


‘3F^ 

3Fi  1 

1 

3(j) 

3(|)2  I 

• • • 

1 3<l)s. 

If  any  of  the  above  sufficient  conditions  are  satisfied, 
then  the  parameter  vector  (|)  is  locally  identifiable.  Note 
that  the  first  condition  is  independent  of  input;  i.e.,  no 
D(({)),  or  u(t)  are  required.  Therefore  condition  one 
allows  checking  for  local  identifiability  independent  of 
input.  The  second  condition  is  independent  of  the  state 
x(t),  but  requires  a nonzero  input  u(t).  This  second  con- 
dition agrees  with  that  derived  recently  by  Grewal  and 
Glover  (Ref  12)  using  linear  system  theory.  The  second  con- 
dition may  be  easier  to  apply  than  condition  one  in  some 
cases,  but  clearly  can  only  be  applied  where  B(<j>)  ^ 0. 

Note  that  these  two  conditions  are  analogous  to  separating 
a linear  system  solution  into  the  zero  input  solution  and 
zero  initial-state  solution,  with  the  total  solution  as  the 
sum. 
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The  above  results  apply  to  either  discrete  or  continu 
ous  observations.  To  test  for  global  identifiability , one 
would  have  to  assure  that  F was  norm-coercive  or  continu- 
ously differentiable  on  an  open  and  path-connected 

Consider  the  following  examples: 


Example  1 

= - ({|)j^  + 4>2)x^(t)  + 4>2  ^^(t)  + u(t) 
x^Ct)  = x^(t)  - (l)^  x^Ct) 

with  observation  process 

y(t)  = Xj^(t) 


This  is  linear  time  invariant 

C 1 0 ] 


A(4>)  = 


C(<1))  = 


’2 
"*3 


system  with 


D(4>)  = [o] 


Therefore  we  can  use  the  results  of  Section  3.4.  First  we 
can  check  if  the  system  is  identifiable  with  zero  input. 
From  Section  3.4  we  use  the  first  sufficient  condition: 
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1 


0 


F^  = 


C(<|)) 

C(4.)A((1)) 

C((|))A^((l.) 

Lc((1))a^(4.)J 


- ((|)i+4i2) 

+ <l>  2^  ^ + 2 4>  3^  <t>  2 + 2 (j>  2 ^ 

3 2 2 

-4)1  -34>1  4)2“54>i4)2 


-4>i4>2+4>2  “<I’2'^3 
24>i^4'2‘*’<J’i<I>2*I>3 

-4'2^4>3+*2‘*’3^ 


94> 


■ 0 
-1 


+24)i+24>2 

-34>i^-64)3^4)2-54)2^ 


0 

0 

■<<>2 

44>i4>2+^2'^3 


0 

-1 

24>i+44>2 

-34)j2_i0^£^^ 

"2<I>2<I’3"24>2^ 


0 

1 


0 

0 

0 


-4>2+2<f>2~^3 
24’i^+4>i4>3~24)24>3+4)3^  -4>2' 


0 

0 


■^2 

*^1*^2”"*^2  '*'2*I’2*^' 


Clearly  this  has  rank  three  for  4*2  ^ 0*  Thus  the  system  is 

O 

locally  identifiable  on  {R  {4>2  “ We  can  also  check 

the  second  condition  developed  in  Section  3.4.  The  second 


r condition  is  to  determine  if  9F2/84)  has  rank  s where 


0 

0 


0 

0 


0 

0 


-1  -1 
2<frj_+2^2  ~ 

-3*^^-6(>j^^$^-5jj>^'^'^-3<|»3^2-10(|)^4.2-2(|)2(|)3-3(l)2^ 


0 

0 


Again  this  has  rank  three  for  ^ Q.  The  above  result 
agrees  with  that  obtained  in  Ref  25. 


Example  2 i 

i 

il(t)  = -<1,2X2  (t)  I 

X2(t)  = -<J,j^Xj^^(t)-<})2X2(t)+u^(t) 

y(t)  = Xj^(t) 

Using  Theorem  3.4,  we  form  F. 

Fo  = x^ (t) 

^1  = -'*'2^2^^^ 

^2  “ “4>2  (“4>iXj^^  (t) -4,2X2  (t)+u^  (t)  ) 

o o 

F3  = +<f,2<|,23Xj^  X3(t)+(|,2  X2  (t)-2<{,2u(t)u(t) 

= -3(|,2^<I>iXj^2  (t)  X2  (t)  +<J)2^  (-4'lXj^^  (t)  -<I>2^2  ^ 

-2<j»2u(t)u  (t) 

= -3(|,3^ <(,2^X3^^  (t)x2(t)-<l,i<t,2^X3^^(t)-<|)2^X2(t) +<1,2  V(t) 

-2(|,2U  (t)  u (t) 


80  that 


x^(t) 

^ -4*2*2 

4»l<j>2Xi^  (t)  +())2^X2  (t)  — <|)2U^  (t) 

-3(|)3^(()2^Xj^^(t)x2(t)-(|)j^<|)2^X3^^(t)-())2^X2(t)+())2^u^t)-2(J)2u(t)u(t) 


3F 

^ = 


0 

0 

-3<|)2^Xj^2  (t)  X2  (t)  “<l>2^^1^  (t) 


0 

-X2(t) 

***1*1^  (t)  +24I2X2  (t)  -u2  (t) 
-6(J)j^(1)2X^^  (t)  X2  (t)  -2(J)j^(j)2X2^^(t) 
-3(p2^X2(t)+2(p2U^(t)-2u(t)u(t) 


This  matrix  has  rank  two  for  (J>2  5^  0.  Therefore  the  param- 

2 

eter  vector  (J>  is  locally  identifiable  for  4)  e $ = {R  ~{<|>2=0^^* 

3.5  Stochastic  Identifiability  of  Parameters  of  Nonlinear 
Dynamical  Systems 

In  this  section  we  are  interested  in  extending  the  pre- 
vious results  to  the  case  in  which  the  observations  are 
corrupted  by  additive  noise.  To  do  this  we  must  first 
define  identifiability  in  terms  of  output  distinguishability 
in  a stochastic  sense. 

The  nonlinear  system  is  modeled  by 

x(t)  = f(t,  x(t),  u(t),  (|))  x(t  ) = X (3.5.1) 

o o 

for  all  t.  The  discrete  observation  process  is  corrupted 
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by  a zero-mean  white  Gaussian  sequence  {vCt^)};  i.e,, 

= y(t^)  + v(t^)  i = 1,  2,  • • • k (3.5.2J, 

where 

y(tj^)  = h(t^,  x(t^)  , 4»)  (3.5.3) 

The  objective  is  to  find  a stochastic  analog  to  the 
sufficient  conditions  for  identifying  (j)  that  were  found  for 
^he  deterministic  case.  Tse  and  Anton  (50)  define  sto- 
chastic identif iability  in  terms  of  consistency  in  prob- 
ability. We  will  use  this  definition  which  is  independent 
of  the  estimation  method  used  for  evaluating  the  parameters. 

Definition;  A sequence  of  estimates  _ ^^2  •••  n ••• 

is  said  to  be  consistent  in  probability  if  for  any  6,  e 
arbitrarily  small,  there  exists  an  N(e,  6)  such  that  for 
n > N (e , 6) 

Un  - <l>ol  I > 6)  < e (3.5.4) 

where  <l>^  is  the  true  value  of  the  parameter.  In  other 
words,  the  sequence  of  estimates  converges  in  probability 
to  the  true  value  <|)q  as  n ->■  ®. 

Definition ; The  parameter  (p^  is  said  to  be  identifiable 
if  there  exists  a sequence  of  estimates  which  is 

consistent  in  probability. 


Definltloni  If  a sequence  of  observations  = (z(t^)}, 

BMde,  then  the  maximum  likelihood 
estimate  of  £ is  ven  by 

♦ - max  p(Z.  I ♦)  (3.5.5) 

With  successive  applications  of  Bayes  Rule,  an  alternate 
expression  for  p(Zj^  I <^)  can  be  found. 

p(Zjj  I ♦)  “ p(z(tj^),  , •••  I 4») 

= P(z(tj^>  I I 

= p(z(tj^)  I ♦)p(z(tj^_j^)  I \_2>  I ^ 

e 

k 

= ir  p(z  (t . ) I Z , ()>) 
j=l  ^ 

(3.5.6) 

The  consistency  of  the  maximum  likelihood  estimator  has  been 
considered  by  several  people.  Cramer  (Ref  7)  and  Wald 
(Ref  53)  showed  it  to  be  consistent  under  independent  obser- 
vations for  reasonably  general  conditions.  Maybeck  (Ref  33) 
and  Tse  and  Anton  (Ref  50)  extend  this  work  to  dependent 
observation  sequences.  The  following  theorem  gives  suffi- 
cient conditions  for  the  maximum  likelihood  estimate  to  be 
consistent  in  probability. 

Theorem  3.9 

A 

Let  <j)(Zj^)  denote  the  maximum  likelihood  estimate  of 
4>  e 4>  under  the  observation  sequence  Z . Then  the  estimate 


<KZj^)  is  consistent  in  probability  (i.e.,  cj)  is  identifiable) 
if  for  all  (|>2  e <l>j^  ^ <J>2/  there  exists  an  infinite 
set  S C (where  I"*"  is  the  set  of  non-negative  integers) 
such  that,  for  all  )c  e S,  the  inequality 

*^1^  ^ P^^k  I ^k-l' 

is  satisfied. 

Theorem  3.9  is  proved  in  Ref  50  under  some  rather 
general  assumptions  concerning  the  joint  density  function 
P(Zj^  I «l»). 

Now  return  to  the  nonlinear  model  described  by 
Eq  (3.4.1)  and  the  discrete  observation  process 

z(t^)  = y(t^)  + v(t^)  i = 1,  2,  • • • k (3.5.7) 

with  v(t^)  an  independent  Gaussian  zero-mean  white  noise 
sequence  with 

E{v(t . ) v(t . ) } = R(t.)  for  t.  = t. 

ID  i 13 

= 0 for  t^  ^ tj  (3.5.8) 

An  observation  vector  Zj^  is  formed  consisting  of  a sequence 
of  observations 

Zk  = |^z‘^(tj^)z'^(t2)  • • • z'^(tj^)]'^  (3.5.9) 

with  t^  e jt^,  tjj  for  all  i. 

z(tj^)  and  v(t^)  are  m-vectors  so  that  Zj^  has  dimension  km. 


The  joint  probability  density  function  of  the  m com- 


ponents of  z (t j ) is  Gaussian  and  is  given  by 


p(z(tj)  |4)) 


T 

exp{-  |Jz(t^)-y(tj)j 

[R(t^)]”^[z(tj)-y(tj)]  } 


(3.5.10) 


Note  that  the  mean  of  p(z(tj)  | <J>)  is  given  by 
E{z(tj)}  = E{y(tj)  + v(tj)}  = y(tj).  We  should  also  note 

that  the  probability  density  of  Eq  (3.5.10)  is  conditioned 
on  y(tj)  which  is  dependent  on  the  model  differential 
equation  and  output  model  with  initial  conditions  and 
input  u.  Also  since  the  observation  noise  is  independent 
in  time 

p(Zj^  I <(i)  = p(z(t^),  2(^2)/  • • * » I <f») 


= -n  p(z(t.)  I (|)) 
i=l 


(3.5.11) 


where  p(z(t^)  | ((>)  is  given  by  Eq  (3.5.10). 

Note  that  an  equivalent'*expression  of  Eq  (3.5.11)  is 


P(Z^  I (|))  = 


, km/2 


«p{-  - y]''  a-1[z^  - y]  ) 


(3.5.12) 


where 


T iT 

z (t,^)] 

T 1 ^ 

y t^k)] 


mk  X mk 

Now  suppose  we  repeat  the  process  on  [t^f  many 

times  so  that  we  can  obtain  independent  repetitions  of  the 
observation  process  Zj^.  We  can  form  a vector  Zj^  from  these 
independent  trials  of  Zj^,  i.e., 

= [“‘i’'  4’’  • • • ’’ 

Note  that  we  have  required  that  the  observation  noise  be 
independent  between  successive  trials  of  Zj^,  so  that 

Ir  Jr 

E{v  (tj)v  (tj)}  = 0 for  trials  i and  i,  i.  ^ Si.  We  can  make 
a maximtim  likelihood  estimate  of  <p  based  on  this  observation 
sequence  as  follows: 

$(N)  = max  P(Z^^,  Z^.^,  • • • Z ” | ()>) 

4>e*  K K K 

N ^ 

= max  V pCZ^"*  I 4>)  (3.5.14) 

♦e*  j=i 

where  p(Zj^^  | ♦)  is  given  by  Eq  (3.5.12). 

From  Theorem  3.9,  we  know  that  a sufficient  condition 
for  ^ to  be  identifiable  is  that  for  all  4>2  c 4>  with 


Zj^  = [^z'^(tj^)z'’^(t2)  • • • 

Y = j^y^(tj^)y^(t2)  • • • 


and 


R(t.)  0 

1 R(t2) 


A » 


R(tj^) 
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^2'  existiL  an  infinite  set  S C l"*"  such  that,  for 

all  j e S,  the  inequality 

I *1*  ^ P^^k^  I ^2^  (3.5.15) 

is  satisfied.  Note  that  since  we  have  independent  trials 
of  the  observation  process  Zj^,  the  following  is  equivalent 
to  the  inequality  Eq  (3.5.15) 

P(Zj5^  I ^i)  p(Zk^  I 4*2^  (3.5.16) 

where  P(Zj^^  | <|))  is  given  by  Eq  (3.5.12). 

Prom  the  form  of  p(Zj^^  | (j>)  in  Eq  (3.5.12),  it  is 
apparent  that  the  inequality  Eq  (3.5.16)  is  equivalent  to 

!<  [Z|^^  - a"^]z|^^  - Y(4'2>]  (3.5.17) 

Consequently  a sufficient  condition  for  Eq  (3.5.17)  to  hold 
is  that 

■ Y(4)j^)  ^ Y((})2) 

for  all  ((>.,,  <1)_  e ♦ with  <j>  ^ <|>  . 

^ X 2 


But  this  is  just  a statement  of  (|)  being  identifiable  on  $ 

under  the  deterministic  observation  process  \ 


I 


The  foregoing  discussion  is  proof  of  the  following 


r 

■i 

I 

[I 

1 

f 

i 

theorem: 

Theorem  3.10 

Given  the  system  modeled  by 

x(t)  = f(t,  x(t),  u(t),  (j))  x(t  ) = X 

o o 

and  the  discrete  observation  process 

z(t.)  = y(t.)  + v(t.)  i = 1,  2,  • • • k 

I XXI 

! 

where 

y (t^)  = h(t^,  x(t^)  , (()) 

and  v(tj^)  is  a zero  mean  white  Gaussian  sequence;  then  the 
parameter  vector  (()  is  identifiable  on  $ (in  the  sense  of  the 

1 

definition  following  Eq  (3.5.4))  using  the  maximum  likelihood  j 

t 
4 

t 

estimate  1 

^ N ^ f 

4>(N)  = max  Ti  p(Z  I <j>) 

({ie4>  j=i  ^ ‘ 

if  and  only  if  (|>  is  identifiable  on  $ under  the  determin- 
istic observation  procese 

Y = [y'^(ti)y'^(t2)  • • • y'*'(tk)] 

Note  that  we  have  required  the  observation  noise 
(vCt^)}  to  be  independent  in  time  and  independent  between 
successive  trials  of  the  observation  process  Zj^.  That  is, 

I 

t 

1 
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for  the  trial  of  the  observation  process;  i.e., 
we  require 

E{v^  (tj^)v^  (tjj)  } = 0 for  i A 

emd  for  the  trials  i and  t of  the  observation  process,  i.e., 
Zj^  and  Zj^  , we  require 

E{v^  (t^)v*' (tj^)  } = 0 for  i j = k,  and  j k. 

Also,  if  $ is  an  open  neighborhood  of  (p,  then  Theorems  3.1, 

3.2,  and  3.4  can  be  applied  to  test  for  deterministic  local 

s 

identif lability.  If  $ C R , Theorems  3.6  and  3.7  can  be 
applied  to  test  for  deterministic  global  identifiability. 

The  above  theorem  can  be  applied  to  continuous  obser- 
vation processes  on  the  interval  t^l*  This  can  be  seen 

by  simply  considering  k > ^ time  points  on  this  interval 
and  then  apply  the  above  development. 

It  is  emphasized  that  Theorem  3.10  can  be  applied  only 
to  cases  of  additive  white  noise  on  the  observation  where 
independent  repetitions  of  the  observation  can  be  obtained. 

If  there  were  driving  process  noise,  the  results  would  be 

I 

t‘  substantially  different.  Generally,  deterministic  iden- 

tifiability is  a necessary  but  not  sufficient  condition  for 
stochastic  identifiability. 
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3.6  Sunmiary 

In  this  chapter,  we  have  developed  sufficient  con- 
ditions for  identifiability  of  nonlinear  systems. 

Local  Identifiability.  Theorems  3.1  and  3.2  state 
sufficient  conditions  for  local  identifiability  under  dis- 
crete observations.  Theorem  3.4  uses  a recursion  relation 
to  present  alternative  sufficient  conditions  under  discrete 
observations.  Theorems  3.3  and  3,5  state  corresponding 
sufficient  conditions  for  local  identifiability  under  con- 


Chapter  IV 


COMPUTATIONAL  TECHNIQUES  FOR  PARAMETER  IDENTIFICATION 

In  the  previous  chapter,  we  defined  local  and  global 
identifiability  in  terms  of  a one-to-one  correspondence 
between  the  parameter  vector  and  the  observed  data.  These 
definitions  are  independent  of  the  computational  technique 
used  for  identification.  Typically  the  parameter  identi- 
fication problem  involves  estimating  certain  parameters 
based  on  experimental  data.  The  general  approach  is  to 
establish  a criterion  or  cost  functional  which  is  used  to 
choose  the  best  values  of  the  parameter  vector  elements. 

We  will  be  concerned  with  parameter  identification  in  which 
the  least  squares  cost  functional  is  applied.  Again  we 
consider  the  class  of  systems  modeled  by 

x(t)  = f(t,  x(t),  u(t),  (j))  = Xq  (4.1) 

for  t e [tQ,  tf]  with  observations 

y(t^)  = h(t^,  x(t^),  (J>)  i = 1,  2,  • • • , k (4.2) 

For  a given  input  function  u,  initial  conditions  x^  and 
observed  experimental  data  z(t^)  for  all  tj^  e *"f^  ' 

we  want  to  choose  (|>  such  that 

T 

Z(t^)  - y(t^)j  l^z(ti) 

is  minimized. 
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In  this  chapter  we  will  develop  some  theorems  con- 
cerning minimization  of  cost  functionals,  relate  identi- 
fiability  with  the  problem  of  minimizing  the  least  squares 
cost  functional,  and  discuss  two  iterative  computational 
techniques  for  determining  a value  of  <{>  which  minimizes 
the  cost  functional. 

4.1  Minimization  of  Cost  Functionals 

In  this  section  we  state  some  definitions  and  theorems 
which  apply  to  the  problem  of  finding  minimum  values  of 
cost  fxinctions.  The  following  definitions  are  taken  from 
Ref  39: 

Definition;  Let  g:DCR^  -^R^.  A point  x*  e D is 
a local  minimizer  of  g if  there  is  an  open  neighborhood  S 
of  X*  such  that  for  all  x e sOd,  x 5^  x* 

g(x)  > g(x*)  (4.1.1) 

X*  is  a proper  local  minimizer  of  g if  strict  inequality 
holds  in  Eq  (4.1.1).  If  Eq  (4.1.1)  holds  for  all  x e C D 
and  X*  e Dq,  then  x*  is  a global  minimizer  of  g on  Dq. 

Definition:  A point  x*  e Int  D (denotes  interior  of  D) 
is  a critical  point  of  g ; D C r”  R^  if  g has  a derivative 
at  X*  and  9g(x*)/9x  (denoted  by  g'(x*))  equals  zero. 

Lemma  4.1 

n X 

If  X*  e Int  D is  a local  minimizer  ofg;DCR  -^R, 

2 2 

and  g' (x*)  and  9 g(x*)/9x  (denoted  by  g''(x*))  exist,  then 
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g' (x*)  = 0 and  g''{x*)  is  positive  semide finite. 

Lennna  4.2 

Let  g : D C r’™  -*■  and  suppose  g'  ' (x)  exists  at 
X*  e Int  D.  Then  x*  is  a proper  local  minimizer  of  g if 

(1)  g'  (X*)  = 0 

and 

(2)  g''(x*)  is  positive  definite. 

Lemma  4.1  and  4.2  are  the  well  known  necessary  and  suffi- 
cient conditions  for  g to  have  a minimum  at  x*. 

Definition;  A functional  g : D C R^  ^ R^  is  convex 
on  a convex  set  C D if  for  all  x,  y,  e and  0 < a < 1 

g(a  X + (l-a)y)  £ a g(x)  + (l-a)g(y) 

The  functional  is  strictly  convex  on  if  strict  inequality 
holds  whenever  x ^ y. 

Theorem  4 . 1 

Assume  that  for  g : D C r"  -»■  R^,  g*  ' (x)  exist  for  each 
point  of  a convex  set  D^^C  D.  Then  g is  convex  on  if  and 
only  if  g' ' (x)  is  positive  semidefinite  for  all  x e D^. 
Moreover,  g is  strictly  convex  on  if  g' ' (x)  is  positive 
definite  for  all  x c D^. 

Proof:  See  Ref  39:87. 
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From  the  Theorem  4.1,  we  see  for  example  that  the  functional 
T 

g(x)  = X Ax  is  convex  if  and  orly  if  A is  positive  semi- 
definite  and  strictly  convex  if  A is  positive  definite. 

Theorem  4 . 2 

n X 

Suppose  that  g ; D C R R is  convex  and  differ- 
entiable on  an  open  convex  set  C D.  Then  x*  e D is  a 
critical  point  of  g if  and  only  if  x*  is  a global  minimizer 
on  Dq.  Moreover,  if  g is  strictly  convex  on  D^,  there  is 
at  most  one  critical  point  in  D^. 

Proof:  See  Ref  39:101. 

From  Theorems  4.1  and  4.2,  we  can  state  the  following 
corollary. 

Corollary  4.1 

Given  g : D C r”  ^ R^.  If  g' ' (x)  is  positive  definite 
at  each  point  of  an  open  convex  set  C D,  then  g has  at 
most  one  local  minimizer  in  D^.  That  is,  if  a critical 
point  exists  in  D^,  it  is  the  only  critical  point  in  and 
is  a global  minimizer. 

Now  consider  a particular  class  of  functionals  which 
we  will  refer  to  as  quadratic.  Phillipson  (42)  has  devel- 
oped some  useful  fundamental  results  pertaining  to  quad- 
ratic functionals  defined  on  a Hilbert  space.  For  our 
purposes  we  can  assume  the  Hilbert  space  to  be  r’^  with 
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inner  product  (x,  y)  = x y and  with  the  Euclidean  norm. 
Let  X,  y,  z e r”  and  consider  the  following  definitions 
given  by  Phillipson; 

S,  (x)  is  a continuous  linear  functional  on  if 
I A (x) I £ kl |x| I k < “ 

and  for  a,  g e R^ 


A (ax  + gy)  = aA(x)  + gA(y) 

q(x,  y)  is  a continuous  bilinear  functional  on  r”  if 


|q(x,  y)  I i k|  |xl  I I |yl  I k < ~ 

and  for  a,  g e R^ 

q(ax  + gy,  z)  = aq(x,  z)  + 3q(y,  z) 
q(x,  ay  + gz)  = aq(x,  y)  + gq(x,  z) 

h(x)  is  a quadratic  functional  on  R*'  with  the  property 


Proof;  From  the  definition  of  a convex  functional,  we  need 
to  show  that 


h(ax  + (l-a)y)  < ah(x)  + (l-a)h(y)  0 < a < 1 
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From  the  definition  of  h(x) 


(l-a)h(y)  + ah(x)  = (l-a)q(y,y)  + aq(x,x) 

- 2[(l-a)il(y)  + ai(x)] 

= q(y»y)  + «q(y/  x-y)  + aq(x-y,  y) 

+ aq(x-y,  x-y)  - 2££(l-a)y  + ax] 

2 

Since  aq(x-y,  x-y)  > o q(x-y,  x-y)  we  obtain 

(l-a)h(y)  + ah(x)  > q(y,y)  + aq(y,  x-y)  + aq(x-y,  y) 

+ a q(x-y,  x-y)  - 2£C(l-a)y  + ax] 
= q(y  + a (x-y),  y + a (x-y)) 

- 2£[ (l-a)y  + ax] 

= h( (l-a)y  + ax) 


Q.E.D. 


Theorem  4.3 


Let  h(x)  be  a quadratic  functional  defined  on  a closed, 
convex  set  D C r”.  Then  there  exists  an  x*  s D which  is  a 
unique  global  minimizer  on  D.  That  is,  there  is  a unique 
X*  such  that 


h(x)  > h(x*) 


for  all  X c D. 


Proof;  The  existence  is  proved  in  Appendix  2.7  of  Ref  42. 
Uniqueness  follows  quickly  from  h(x)  being  strictly  convex. 
Suppose  we  assume  there  is  an  x**  ^ x*  with  x**  e D such 
that 


h(x**)  = h(x*)  “ h^.„ 


Prom  Lemma  4.3,  we  know 


h(ax*  + (l-a)x**)  < ah(x*)  + (l-a)h(x**)  = h . 

mi 

Now  let  ax*  + (l-a)x**  = y e D. 

Then 


h(y) 


which  is  a contradiction. 

Q.E.D. 

Theorem  4.3  is  a fundamental  result  which  assures  us 
that  a unique  global  minimizer  of  h(x)  exists  on  a closed, 
convex  set  D C R*^.  Two  additional  characteristics  of  the 
global  minimizer  of  h (x)  are  contained  in  the  following 
theorems . 

Theorem  4.4 

The  unique  global  minimizer  x*  for  which 
h(x*)  ^ h(x)  for  all  x e D 

is  characterized  by 

q(x*,  x-x*)  ^ Jl(x-x*) 

for  all  X e D. 

Proof:  See  Ref  42:20. 
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Theorem  4.5 


For  D = r”,  then  the  unique  global  minimizer,  x*,  is 
characterized  by 

q(x*,  x)  - A(x)  =0  for  all  x e r”. 

Proof:  See  Ref  42:21. 

In  the  next  section  we  will  apply  these  results  to  the 
least  squares  cost  functional. 

4,2  Identif iability  and  Minimization  of  Cost  Functionals 

In  this  section,  v/e  want  to  relate  identif  iability  with 
the  problem  of  minimizing  the  least  squares  cost  functional 
given  by  Eq  (4.3).  Recall  that  the  discrete  observation 
process  may  be  represented  by  a vector 

T 

Y(«j))  = Cy'^(tj^)y'^(t2)  • • • y^(tj,)]  (4,2.1) 

and  the  corresponding  experimental  data  by 

T 

Z = Cz'^(tj^)z'^(t2)  • • • z^(tj^)]  (4.2.2) 

Then  the  least  squares  cost  functional  takes  the  form 

J(<J>)  = [z  - Y((f.)]'^  [z  - Y((|))]  (4.2.3) 

We  also  note  that  J(4>)  is  a quadratic  functional  as  previ- 
ously defined  when  viewed  as  a function  of  Y,  That  is 


J : n 


C R -*■ 


is  a quadratic  functional  since 


J(Y)  = y'^Y  - 2[z'^Y  - i z'^’z] 
= q(Y,  Y)  - 2Z(Y) 


Theorem  4.6 

Let  the  class  of  models  be  systems  described  by 

x(t)  = f(t,  x(t),  u(t),  ({))  t e [t^,  tj] 

= '‘o 

and  the  discrete  observation  process  be  described  by  Y(((>) 
with  the  corresponding  experimental  data  Z.  Further  let 
Y{(|))  ; $ C n C R™^  be  a mapping  such  that  fi  is  a 

closed  convex  set.  Then  the  least  squares  cost  functional 
J(<|))  has  a unique  global  minimizer,  4>*e$  if  the  system  is 
globally  identifiable  on 

Proof:  Since  J (Y)  is  a quadratic  functional,  we  know  from 

Theorem  4.3  that  there  is  a unique  Y*  e which  is  a global 
minimizer  of  J on  fi.  Also  if  the  system  is  globally 
identifiable  on  then«» 

Y(4>i)  Y(4>2) 

for  all  ({i^,  4>2f  e with  ^ <p2»  Therefore  there  is  a 
unique  which  maps  to  Y*. 

Q.E.D. 

Note  that  the  requirement  of  a closed  set  is  needed 
for  existence.  The  uniqueness  follows  for  any  convex  set  n 
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from  the  fact  that  J is  convex  on  n.  Although  we  proved 

convexity  for  any  quadratic  functional  in  Lemma  4.3,  we 

2 2 

could  also  show  J (Y)  is  convex  by  observing  that  d J(Y)/dY 
is  positive  definite  on  Si  and  from  Theorem  4.1  this  implies 
convexity.  Also,  the  uniqueness  of  a global  minimizer,  if 
it  exists,  on  an  open  convex  set  follows  immediately  from 
Corollary  4.1. 

Theorem  4 . 7 

If  (j»  is  mapped  onto  by  Y,  that  is  Y ; <I>  C R®  -»■  R."*^, 
then  the  unique  global  minimizer  (fi*  of  J is  characterized  by 

Y(4.*)  = Z 

Proof;  From  Theorem  4.5  we  know  that 

Y*^  (<!,*)  Y(,j,)  - 2[z^Y(4,)  - 7 Z^Z]  = 0 
mV 

for  all  Y e R . At  Y(<j>)  = Y((|)*),  this  implies 


Y(<|)*)  = Z 

Q.E.D. 

The  above  theorem  can  be  argued  intuitively,  since  if 
one  has  the  choice  of  any  Y(<J))  e R*”^,  it  is  clear  that  the 
least  squares  cost  functional  will  have  its  minimum  value 
of  zero  at  Y(<j)*)  = Z. 

The  following  theorem  relates  local  identif lability  to 


local  minimum  of  J(4>) 


Theorem  4,8 


If  e 4 is  a critical  point  of  the  least  squares 
cost  functional  J((|))  and  3Y((p*)/3(p  (denoted  by  Y*((p*))  has 
rank  s,  then  (j)*  is  a unique  local  minimizer  on  an  open 
neighborhood  S of  (J>*. 

Proof:  If  <p*  is  a critical  point,  then 

3J(<|)*)  .T 
0 

or  from  Eq  (4.2.3) 

T 

“2[y'  ((J)*)]  [Z  - Y(<|)*)]  = 0 (4.2.4) 

Now  since  Y'  ((j)*)  has  rank  s,  there  exists  a nonsingular 
mk  X mk  matrix  Q such  that  Yj^'  ((J>*)  given  by 

Yl'  (<!»*)  = Q Y'  (({.*) 

is  in  Hermite  normal  form  (see  Ref  37:54).  The  form  of 
Yj^'  (4»*)  is  uniquely  determined  by  Y'  (((>*)  and  will  be 
characterized  as  follows: 

Yl'  (<!>*)  = 

where  Y2' (<!>*)  is  an  s x s nonsingular  matrix. 

Eq  (4.2.4)  is  equivalent  to  the  following 
T T T ""X 

-2LY'((|.*)]  qCq]  [z  - Y (<!)*)]  »0  (4.2.5) 

i 

where  Q is  the  nonsingular  mk  x mk  matrix  yielding  the 
Hermite  normal  form  of  Y*  (({)*). 
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Define  the  following 


[^3*  (♦*)]’’  = [y*  ((|)*)]^  q’’  = [q  Y'  ((J)*)]'' 
= [Yi’(<|)*)]^ 


Z3  - =[q'^"^[z  - Y(<j)*)] 


Hote  that 


Y,'((|)*)  = 


0 


where  Y^'  (<|)*)  is  s x s and  nonsingular.  Then  Eq  (4.2.5) 


takes  the  form 


rY4'((|>*)  ^ . 

.2  ^ [Zo  - Y,  (())*)]  = 0 

0 ^ ^ 


■Y4*  [Z4  - Y^Cd.*)' 

0 J U5  - ^5(4.*)] 


This  yeilds 


-^Cy^’  (4)*)]^  CZ4  - Y4(4>*)] 


"2[oy  [Zc  - Yc(4>*)l  = 0 (4.2.6) 


The  first  term  of  Eq  (4.2,6)  can  only  have  the  solution 


?4  = Y4(4.*) 


9inoe  Y4'  (4>*)  is  nonsingular. 


• • y "r- 


;•  iv  '{I 


Recall  from  Theorem  3.1  that  if  Y has  a strong  derivative 
at  (|i*  and  Y'  (<|>*)  has  rank  s,  then  the  system  is  locally 
identifiable  at  (j>*.  This  means  Y{()))  is  one-to-one  on  some 


open  neighborhood  of  (J>*  and  therefore  only  (p*  e Sj^ 
satisfies  Y^(4>*)  = Z^. 

2 2 

We  will  use  the  above  result  to  show  that  3 J/d<p 

is  positive  definite  on  some  neighborhood  of  <>*.  We  obtain 

2 2 

an  expression  for  3 J (({>*) /3<()  expressed  in  the  coordinate 
system  resulting  from  the  previous  transformation  by 
taking  the  derivative  of  Eq  (4.2.6);  i.e., 

~ - -2 ( [i/ (♦*)]’■  [z, 

- [Y4' (<))*)]'^  [y4'  (<()*)]}  (4.2.7) 

From  the  previous  discussion  we  know  that  a neighborhood 
S2  = {<}>  < 5)  exists  such  that  for  (j)  e S2 

II  [Y4"  [Z4  - I i < e . 

The  second  term  of  Eq  (4.2.7)  can  be  written  as 

M4((|)*)  = Cy^' (<j>*)]'^  [ Y^'  (<!>*)] 


r8Y4(ti)l 

T 

■ 3Y4(tj^)  ■ 

3<}> 

3(j) 
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(<(>*)  is  a real  synunetric  matrix.  We  proved  in  Theorem  2.3 


that  is  nonsingular  if  and  only  if  ((j)*)  has 

rank  s.  Therefore  (<(>*)  is  nonsingular.  Also,  a real 

symmetric  matrix  of  the  form  M^((l))  = A A is  positive 

definite  if  it  is  nonsingular  (see,  e.g..  Ref  20:38). 

Therefore  there  is  an  open  neighborhood  S of  (J>*  where 

2 2 

M^(^)  is  positive  definite  and  such  that  3 J((f))/3(|) 
is  positive  definite  for  <t>  e S. 

Then  from  Corollary  4.1  we  know  that  4>*  is  a unique 
local  minimizer  of  J on  S. 

Q.E.D. 

Theorem  4.9 

If  (|)*  is  a critical  point  of  the  least  squares  cost 
functional  and  Y' (({)*)  has  rank  s for  a system  which 

is  globally  identifiable  on  an  open  set  $ C D,  then  (f>*  is 
a unique  global  minimizer  of  J on  $. 


Proof;  Following  the  same  procedure  as  the  proof  for 
Theorem  4.8,  we  can  show  that  4>*  is  the  only  critical  point 
of  J on  $.  That  is,  the  necessary  condition  for  a critical 
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point,  9J(<|)*) /3<|>  = 0,  has  a unique  solution  Since 

the  system  is  globally  identifiable  on  4>,  then  there  is 
only  one  (p*  which  can  be  mapped  to  Y(4>*).  To  show  that  41* 
is  a minimizer  of  J,  we  can  use  exactly  the  same  arguments 
as  in  the  previous  theorem  to  show  that  4>*  is  a local  mini- 
mizer on  some  open  neighborhood  of  4)*.  Then  since  4>*  is 
the  only  critical  point  on  it  necessarily  must  be  a 
global  minimizer  on 

Q.E.D. 

The  above  discussion  illustrates  that  identifiability 
and  minimization  of  the  least  squares  cost  functional  are 
intimately  related.  Theorem  4.8  has  given  us  a very  useful 
result.  If  through  some  manner  we  can  find  a critical 
point  4>*  oi  the  least  squares  cost  functional  J,  then  if 
Y'  (4>*)  has  rank  s,  we  not  only  assure  local  identifiability 
but  also  assure  that  4>*  is  a local  minimizer  of  J.  We 
also  know  that  Y'  (4>)  has  rank  s if  and  only  if 

M(4>)  = [y*  (4-*)]'^  [y*  (4)*)] 

is  nonsingular. 

In  the  next  section  we  will  discuss  computation 

T 

techniques  for  solving  the  equation  9J/94>  = 0 , the  solution 
of  which,  of  course,  defines  a critical  point. 

4.3  Computational  Techniques 

Two  computational  techniques  will  be  briefly  presented 
for  minimizing  the  least  squares  cost  functional.  A basic 
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gradient  method  will  be  developed,  followed  by  a quasi- 
linearization technique. 


4.3.1  Gradient  Method 

The  gradient  method  is  an  iterative  technique  in  which 
each  step  is  in  the  direction  of  steepest  decent  toward  the 
minimum  of  the  cost  functional.  We  know  from  Lemma  4.1  that 
a necessary  condition  for  a minimum  is 


_ ^ 
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make  the  largest  step  towards  the  minimum  of  J.  The  step 
size  is  constrained  by 

,T 


Cm'^]  [ 


(4.3.S) 


To  exploit  the  unconstrained  optimization  technique,  we 
will  want  to  minimize 


^ [44'']  - q2) 


(4.3.6) 


where  v is  a Lagrange  multiplier.  A necessary  condition 
for  a minimum  is  the  first  variation  must  be  zero.  Taking 
the  first  variation  of  Eq  (4.3.6)  and  setting  it  equal  to 
zero  results  in 


(4.3.7) 


From  Eq  (4.3.7),  it  is  apparent  that  we  should  choose 


(4.3 .8) 


.k  . 


where  K is  a positive  scalar. 

Eq  (4.3.8)  tells  us  the  largest  possible  change  in  J 
is  obtained  by  stepping  in  the  direction  of  the  local 
gradient  vector.  Thus  the  iterative  calculation  proceeds 
in  a straightforward  manner.  An  initial  value  of  <p  is 
assumed;  i.e.,  We  calculate  dJ (4il ) /d(J)l  and  then 


determine  <|)  from  Eq  (4.3.8).  That  is 


= 9 +A(J>  - 9 - K 


(4.3.9) 


1; 
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The  calculation  is  repeated  to  converge  on  a critical  point 
of  JC*) . 


r . 


In  some  cases  it  may  be  useful  to  choose  the  (k  + l)-st 
iterate  along  the  direction  of  steepest  descent  in  the  norm 
defined  by 


lIxIL  = (x’'  x)  = 


,k 


where  N is  symmetric  and  positive  definite.  In  this  case 
it  can  be  shown  that  the  direction  of  steepest  descent  of  J 
is  given  by  (see  Ref  20:245) 

rdJ(^1 


- [- 


d(|)* 


Then 


.k+1  _ .k 

tp  - (j) 


-lldJCd)^) 


d(J)^ 


(4.3.10) 


The  gradient  method  is  part  of  a broad  class  of  mini- 
mization algorithms  in  which  the  cost  functional  is 
decreased  at  each  stage;  that  is 


J(<|)^'^^)  < J((j)’^) 


k = 0,  1, 


The  following  lemma  relates  to  the  convergence  of 
such  techniques. 


Lemma  4 . 4 


s 1 

Suppose  J ; D C R R is  differentiable  at 


X e Int(D)  and  that,  for  some  p e R‘ 


' 3x 


p > 0.  Then 
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are  relatively  small;  the  main  computational  problem  is 
the  evaluation  of  the  gradient  vector  at  each  iteration. 
This  technique  has  the  disadvantage  that  it  converges 
slowly  as  the  solution  is  approached.  This  is  because  the 
gradient  becomes  small  as  the  critical  point  is  approached. 
Another  drawback  is  that  the  gradient  technique  is  rela- 
tively inefficient  when  ridges  and  ravines  are  encountered 
in  stepping  toward  the  optimal  solution.  The  technique 
tends  to  "walk"  between  the  sides  of  the  ravine  while 
progressing  slowly  out  of  it. 

4.3.2  Quasilinearization  Method 

The  preceding  section  developed  a direct  computational 
technique.  In  this  section  we  will  briefly  examine  an 
indirect  method  known  as  quasilinearization.  Indirect 
methods  result  from  attacking  the  solution  to  the  two  point 
boundary  problem  which  arise  from  applying  optimization 
theory.  The  quasilinearization  method  is  also  an  iterative 
procedure,  initially  developed  by  Bellman  and  Kalaba. 

Again  we  are  trying  to  minimize  the  least  squares  cost 
functional  — 

N ™ 

J(«f)  = I [z(t.)  - y(t.)J  [z(t.)  - y(t.)]  (4.3.11) 

i=l  ^ ^ ^ ^ 

Applying  the  necessary  condition  for  a critical  point,  we 
set  3J/3<j»  = 0 to  obtain 

N ,n 

I Cy'(t.)]  [z(t.)  - y(t.)]  =0  (4.3.12) 

j j- 


- — ii— 1 'T  rTKitfl 


Asstune  that  we  have  an  initial  approximation  for  the  pareim- 
eter  vector,  which  again  is  denoted  by  (j)  . Then  expand 
<l>)  in  a Taylor  series  about 

y(t^)  = y(t^f  (p^)  +[y'(tj_,  + 


• • • 


(4. 3.13) 

Neglecting  higher  order  terms  and  substituting  Eq  (4.3.13) 
into  Eq  (4.3.12)  yields 
N V T 

I [y' (^if  {z(t^)  - y(ti,  (P^) 

i»l  ^ 

- y*  (t^,  (p  )[(P  - kJ}  = 0 


(4.  3.14) 


or 


I [y(ti,  <p^)3  [y'(ti,  4.’')]  Ll>  - <P^J 

i*=l 

N T 


I [y'(t..  Ip  )]  [z(t.)  - y(t.,  ())^ 

i“l 


)] 


(4.  3.15) 


Assuming  that  the  inverse  of 


N 


M((^^)  = Z Cy'(ti,  4.’')]  [y'(ti,  <p^)J 

i=l 


(4.  3.16) 


exists,  we  can  obtain  from  Eq  (4.3.15) 

♦ - «=  [m((|)’')]"  ^ [y'(t.  ,T)]  [z(t.) 

i«=i  ^ ^ 


y(tj^,  <p^)J 

(4.  3.17) 
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Eq  (4.3.17)  provides  the  basis  for  the  quasilinearization 


iterative  technique  for  converging  to  a critical  point  of 

k k+1 

J.  We  assume  a (|)  , and  then  calculate  (f>  from 

+ [m I [y  (t. , 4,'')]’’ 

i=l  ^ 

[z(t^)  - y(t^,  d>’')] 

(4. 3.18) 

The  quasilinearization  algorithm  often  requires  a very 
good  initial  estimate  of  the  unknown  vector  (}>  in  order  to 
converge  (Ref  44:152).  One  approach  is  to  use  the  gradient 
method  to  obtain  an  estimate  (|)  near  the  solution.  As 
discussed  before,  the  gradient  technique  is  capable  of  con- 
verging to  the  vicinity  of  the  solution  starting  from  a 
poor  estimate.  Once  a good  estimate  is  obtained,  conver- 
gence to  the  solution  can  be  quite  rapid  using  the  quasi- 
linearization method.  On  the  other  hand,  examples  have 
been  given  to  demonstrate  that  the  algorithm  may,  however, 
diverge  when  started  arbitrarily  close  to  a solution  (see 
Ref  4)  . 

It  was  noted  above  that  the  matrix  M((j)  ) must  have  an 
inverse  for  the  quasilinearization  algorithm  to  function. 
This  matrix  may  be  singular,  or  nearly  so,  especially  in 
the  initial  iterations.  As  discussed  earlier,  the  non- 
singularity of  M((|))  is  required  to  assure  local  identifi- 
ability.  We  also  found  that  for  the  least  squares  cost 
functional,  that  the  nonsingularity  of  M(c|))  will  be  suffi- 
cient for  the  solution  to  be  a local  minimizer.  Thus  if 
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the  quasilinearization  technique  converges  to  a solution  of 
3J/3<j)  = 0,  then  both  the  necessary  and  sufficient  conditions 
have  been  satisfied  for  a local  minimum.  Banks  and  Groome 
(Ref  43)  arrive  at  this  conclusion  in  their  study  of  the 
cpnvergence  of  this  algorithm.  We  can  make  this  statement 
based  on  Theorem  4 . 8 which  depends  on  the  least  squares 
cost  functional  and  system  identif lability , but  is  inde- 
pendent of  the  computational  technique  used  to  find  a 
critical  point. 

4.4  Summary 

In  this  chapter  we  have  developed  some  theorems  con- 
cerning the  minimization  of  cost  functions,  with  particular 
emphasis  on  the  least  squared  cost  functional  as  given  in 
^q  (4.3).  In  Theorem  4.8,  it  was  shown  that  if  4>*  c 4>  is  a 
Critical  point  of  the  least  squares  cost  functional,  J(4>)# 
which  is  a function  of  the  observation  process  Y((j>),  then 
is  a unique  local  minimizer  of  J((J>)  if  Y’ (<l>*)  has  rank  s. 
Recalling  the  results  of  Chapter  III,  we  know  that  (J>,  the 

* 

vector  of  s parameters,  is  also  locally  identifiable  at  4) 
if  Y'  (4)*)  has  rank  s.  flieorem  4.8  was  extended  to  apply  to 
fin  open  set  $ by  adding  the  requirement  that  4>  he  globally 
identifiable  on  $ (Theorem  4.9).  Two  computational  tech- 
niques were  presented  as  methods  for  minimizing  the  least 
squares  cost  function.  The  gradient  and  quasilinearization 
teohniques  were  developed  together  with  some  of  the  features 
which  characterize  these  two  computational  procedures. 
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Chapter  V 


OPTIMAL  CONTROL  MODEL  FOR  HUMAN  PERFORMANCE 


The  previous  chapters  developed  sufficient  conditions 
for  establishing  identifiability  of  nonlinear  systems.  1 

Also,  the  use  of  the  least  squares  cost  functional  and  the 

maximum  likelihood  technique  to  estimate  parameters  based  : 

on  experimental  data  was  discussed  in  the  context  of  system 
identifiability.  We  want  to  apply  this  theory  to  the  param- 
eter identification  problem  associated  with  the  optimal  | 

control  model  for  human  performance.  This  model  has  proven 
to  be  an  effective  tool  for  modeling  human  performance  in 
various  tasks  (Refs  5,  20,  21,  22,  23).  Previous  attempts 
to  examine  the  identifiability  of  the  model  parameters  have 
been  limited  to  investigations  using  linear  theory 

(Refs  40,  41).  The  purpose  of  this  chapter  is  to  describe  i 

the  optimal  control  model  for  human  performance  and  place 
the  problem  of  parameter  identifiability  in  the  context  of 
the  nonlinear  theory  previously  presented. 

As  will  be  seen  in*^he  subsequent  development,  the 
basic  structure  of  the  optimal  control  model  results  in  a 
set  of  model  parameters,  the  values  of  which  require  iden- 
tification using  experimental  data.  In  Chapter  VII,  we 
will  use  the  theory  of  Chapter  III  to  investigate  the  iden- 
tifiability of  the  uncertain  model  parameters.  We  will 
also  attempt  to  identify  the  values  of  the  parameters  using 
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the  gradient  technique  discussed  in  Chapter  IV  in  con- 
junction with  experimental  data  obtained  from  a simulator 
system.  The  next  two  chapters  will  develop  the  optimal 
control  model  and  relate  the  model  to  the  simulator  system 
used  for  obtaining  experimental  data. 

5.1  Structure  of  the  Optimal  Control  Model 

The  optimal  control  model  for  human  performance 
assumes  the  operator  performance  is  dictated  by  the  desire 
to  behave  optimally  with  respect  to  a chosen  cost  functional. 
Factors  are  incorporated  into  the  model  to  account  for  the 
human's  inherent  limitations  such  as  time  delay  and  ran- 
domness . 

The  basic  structure  of  the  optimal  control  model  is 
shown  in  Fig.  1.1.  It  is  assumed  that  the  system  dynamics 
can  be  described  by  linear  equations  and  these  dynamics  are 
under  the  influence  of  an  operator  control  input  as  well  as 
random  disturbances.  The  operator  observes  a set  of  out- 
puts that  are  a linear  function  of  the  system  state  and 
operator  control;  however,  it  is  assumed  that  the  operator 
perceives  a delayed  noisy  version  of  the  observed  variables. 
The  operator's  control  is  derived  assuming  that  a quadratic 
cost  functional  is  minimized  on  the  basis  of  the  perceived 
information.  The  foregoing  scenario  for  human  performance 
gives  rise  to  the  Kalman  filter,  predictor,  and  optimal  con- 
trol gains  appearing  in  the  model  structure  as  shown  in 
Fig.  I.l. 
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We  will  be  concerned  with  using  the  optimal  control 
model  to  represent  human  performance  in  a tracking  situ- 
ation. That  is,  the  operator  is  attempting  to  accurately 
track  a moving  target  with  the  aid  of  certain  displayed 
information.  In  the  remainder  of  this  chapter,  the  optimal 
control  model  will  be  described  in  some  detail,  with  model 
parameters  being  defined  as  they  arise  in  the  mathematical 
development. 

5.2  Mathematical  Description  of  the  Model 

We  assume  the  system  dynamics  can  be  represented  by  a ! 

stochastic  linear  differential  equation  ; 

j 

dx(t)  = A(t)x(t)dt  + B(t)u(t)dt  + r(t)d3(t)  (5.2.1) 

with  x(t^)  = Xq  and  where 

x(t)  = system  state 

1 

A(t)  = homogeneous  system  dynamics  matrix 
B(t)  = control  input  matrix 

3 

u(t)  = human's  scalar  control  input  ■ 

i 

r(t)  = process  noise  distribution  matrix 

3(t)  = Brownian  motion  of  diffusion  strength  W(t) 
for  all  t 

Proceeding  formally.  Eg  (5.2.1)  can  be  expressed  as 

x(t)  = A(t)x(t)  + B(t)u(t)  + r(t)w(t)  (5.2,2) 

with  x(t^)  = Xq  and 
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where  w(t)  is  zero-mean  white  Gaussian  noise  with 
covariance  kernel 

E{w(tj^)w^(t2)  } = W(tj^)6(tj^  - t2) 

The  state  vector,  x(t),  can  be  partitioned  into  two  parts; 

i.e. , 


X (t)  = [Xj^ 


1 ^£+1 


(5.2.3) 


where  x, , • • • x.  are  the  states  associated  with  the  dis- 
turbance  noise  model  and  are  the  states  of 

the  system  dynamics.  Then  the  matrices  A(t) , F (t) , and 
B(t)  can  be  viewed  as  having  the  following  partitioned 
form 


A(t)  = 


F(t)  = 


B(t)  = 


£x£ 

\ ' 
n 1 

0 

1 

A 1 

A, 

B 1 

d 

nxn 

Fn(t)  - 

£xl 

0 

nxl 

0 ■ 

£xl 

B (t) 
c 

nxl 

(5.2.4) 


The  m variables  displayed  to  the  operator,  y^(t),  are 
linearly  related  to  the  system  state  and  control  as  follows 

y . (t)  = C(t)x(t)  + D(t)u(t)  (5.2.5) 

a 


It  is  assiamed  that  if  a variable  is  displayed  visually  to 
the  operator,  then  rate  of  change  of  that  quantity  is  also 
perceived  by  him  (32) . Therefore  if  a visual  display  is 
utilized,  y^j(t)  includes  all  variables  displayed  explicitly 
plus  the  first  derivatives  of  those  variables. 

To  account  for  the  human  operator's  inherent  limi- 
tations, it  is  assumed  that  the  operator  perceives  a delayed 
and  noisy  version  of  Y^jCt)  given  by 

y(t)  = y^j(t  - T)  + v(t)  (5.2.6) 

where 


T = equivalent  perceptual  delay 
v(t)  = equivalent  observation  noise 


with  v(t)  a zero-mean  white  Gaussian  noise  with  independent 
components  and  variance  kernels 


E{v. (t  )v. (t  ) } = V. (t  )6(t,  - t,) 

XJ.  ^ ^ 


i = 1,  2, 


• • , m 

(5.2.7) 


The  randomness  that  the  human  injects  into  the  system 
(referred  to  as  remnant^,  is  included  in  the  equivalent  obser- 
vation noise,  v(t).  In  other  words,  sources  of  remnant  are 
referred  back  to  the  input  of  the  human  model  and  appear  as 
noise  associated  with  the  observation  of  displayed  vari- 
ables. It  is  very  difficult  to  differentiate  among  various 
remnant  sources,  and  combining  them  into  an  equivalent  obser- 
vation noise  has  been  found  to  be  a valid  procedure  (5) . 
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When  directly  viewing  (t)  , the  magnitude  of  the 
associated  observation  noise  variance  has  been  found 
(Refs  5,  23)  empirically  to  scale  with  the  variance  of 
y^  (t) ; that  is 

V^(t)  = E{y^^(t)}  (5.2.8) 

where  P2»  * * • Pj„  are  constants.  A typical  nominal 

value  for  p^  is  .01  tt  (Ref  23). 

The  perceptual  delay,  t,  is  intended  to  represent  the 
combined  effects  of  sensory  and  information  processing 
delays  within  the  human.  Typically,  t has  a nominal  value 
between  .15  and  .3  seconds  (Refs  5,  23). 

It  is  assumed  that  the  operator  is  attempting  to 
minimize 

1 ^f  T »2 

J (u)  = lira  E{  Q x(t)  + g(t)u  (t)]dt) 

^ (5.2.9) 

based  on  the  perceived  information  y(t).  This  operator 
cost  function  was  chosen  for  its  physical  appeal  when 
applied  to  a tracking  task  and  because  it  leads  to  mathe- 
matical tractability. 

The  weighting  matrix  Q has  the  form 


0 = 


I 

I 


0 


(5.2.10) 


with  this  form  the  first  term  of  the  cost  functional  repre- 


sents a mean  squared  error  criterion  with  specifying  a 
constant  cost  weighting  of  the  system  dynamic  states.  The 
term  g(t)u^(t)  is  used  to  account  for  the  operator's  limi- 
tation on  the  rate  of  control  motion  and,  as  such,  intro- 
duces a "neuro-motor"  lag  in  the  model.  This  will  be  dis- 
cussed further  in  association  with  the  calculation  of  the 

optimal  control  gains.  Note  that  the  cost  functional  could 
. 2 

include  terms  such  as  u (t) . The  form  of  the  cost  func- 
tional is  task-dependent  and  needs  to  be  analyzed  from  the 
standpoint  of  operator  intent  and  system  limitations.  As 
will  be  seen  in  Chapter  VII,  the  experimental  simulator 
used  in  this  research,  as  well  as  the  intent  of  the  oper- 
ators, , is  best  represented  by  a cost  functional  of  the 
form  in  Eq  (5.2.9).  The  assumption  of  t^-*-“>  is  valid  so 
long  as  the  total  tracking  time  t^  is  much  greater  than  the 
system  time  constants. 

The  remaining  elements  of  the  optimal  control  model 
are  the  result  of  solving  a stochastic  regulator  problem 
with  a time  delay  constraint.  It  is  shown  in  Ref  26  that 
the  optimal  control  is  generated  by  a linear  feedback 
system  which  consists  of  a cascade  combination  of  a Kalman 
filter,  a least  mean-squared  predictor,  and  a deterministic 
optimal  controller  gain.  The  Kalman  estimator  generates  a 

A 

best  estimate,  x(t  - x) , from  the  perceived  data.  Then  the 

A . 

predictor  obtains  a least  mean-squared  prediction  x(t  | a) 
from  the  output  of  the  Kalman  filter. 
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The  controller  portion  of  the  model  consists  of  a set 
of  gains  operating  on  the  estimates  of  the  system  states, 
followed  by  a simple  first-order  lag  to  model  neuro- 
muscular delay.  As  shown  in  Ref  26,  the  optimal  gains  are 
those  obtained  from  the  noise— free  optimal  regulator  prob- 
lem and  are  independent  of  the  time  delay,  x.  To  find  the 
feedback  gains,  an  augmented  state  vector  is  formed  as 
follows 


so  that 


Xp(t) 


x(t) 

u(t) 


(5.2.11) 


XR(t)  = A^(t)Xj^(t)  + Bj^u(t) 


= X. 


RO 


(5.2.12) 


where 


A(t) 

I B(t)  “ 

A„(t)  = 

1 

1 

0 

! 0 . 

and 

0 

Br  = 

1 

The  subscript  R is  used  to  denote  the  augmented  state 
equation  which  leads  to  the  Riccati  equation  for  determining 
the  optimal  control.  This  distinguishes  this  augmented 
state  equation  from  one  which  appears  below  with  subscript  a. 
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Prom  the  solution  to  the  linear  regulator  problem, 


the  feedback  gains,  X (t) , are  given  by 

»"=> 


(5.2.13) 


where 

u(t)  = - X(t)  $j^(t)  (5.2.14) 


and  Pp(t)  is  the  solution  co  the  Riccati  equation 


-Pp(t)  = Pp(t)Ap(t)  + Ap(t)Pj^i/..  - Pj^(t)Bj^  -i—  B^P^(t)  + 


g (t) 


(5.2.15) 


where 


Qr  = 


0 

0 


In  the  general  time  varying  case,  the  solution  to 
Eq  (5.2.15)  is  obtained  by  computing  backwards  in  time  from 
t£  to  obtain  Pp(t).  This  is  unrealistic  in  the  human  mod- 
eling context  in  which  the  operator  does  not  know  a priori 
the  time  variation  of  Ap(t).  Therefore  it  will  be  assumed 
that  the  operator  determines  the  gains  adaptively  on  line 
by  using  his  information  at  time  t (i.e.,  Aj^(t)  and  g(t)) 
to  compute  Pj^(t)  backwards  from  t = t^.  That  is,  A^(t)  and 
g(t)  will  be  assigned  their  current  values  and  Pj^(t)  com- 
puted assuming  they  are  constant  on  [t,  t^].  We  assume  t^ 
is  very  large  compared  to  the  system  time  constants,  so 
that  t^-«»  and  therefore  Pj^(t)  is  obtained  by  solving  the 
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steady  state  Riccati  equation  (i.e.,  solve  Eq  (5.2.15)  with 


Pj^(t)  =0).  With  this  adaptive  procedure,  the  solution 
P (t)  of  the  steady  state  Riccati  equation  will  change  at 
each  point  in  time  since  A (t)  and  g(t)  vary  with  time. 

This  will  result  in  the  control  gains,  X(t),  varying  with 
time  (see  Eq  (5.2.13)). 

• 2 

It  was  previously  stated  that  the  term  g(t)u  (t)  in 
the  operator's  cost  functional  results  in  a first  order  lag 
in  the  model.  To  show  how  this  evolves,  it  is  necessary  to 
define  three  new  quantities;  (1)  a first-order  lag  time 
constant,  t ; (2)  an  operator  commanded  control,  u (t) ; and 

Ii  o 

(3)  the  gains  Xj,(t)  which  are  applied  to  the  estimated 

A 

states  x(t)  to  obtain  u (t) . That  is 

c 

u^(t)  = - X (t)x(t)  (5.2.16a) 

c c 

and 

u (t)  = T u(t)  + u(t)  (5.2.16b) 

C 11 

Combining  the  above  two  equations  results  in 

u(t)  = - C x„(t)x(t)  + u(t)]  (5.2.17) 

To  relate  X to  X(t),  we  use  Eq  (5.2.14)  to  show 
c 

fx(t)]  1 

- X(t)|u(t)J  = - [ X^(t)x(t)  - u(t)]  (5.2.18) 

which  implies 

Xi(t)  X^.  (t)  i = 1,  2,  • • • , 1 + n (5.2.19) 

T n c-*- 
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and  the  (Jl+n+l)  component  of  X is  given  by 


With  the  above  definitions,  may  be  given  the  physical 
interpretation  of  a "neuro-muscular"  delay  since  these 
results  yield  a configuration  as  depicted  below: 


x(t) 


-Xc(t) 

1 

o 

p 

1 

U(t) 

T„S+1 

Prom  a physical  standpoint,  this  gives  motivation  to  struc- 
ture the  model  so  that  is  held  constant  (a  nominal  value 
of  .1  second  has  been  found  to  give  good  results 
(Refs  5,  23)).  Note  that  from  (5.2.13)  and  (5.2.20)  one 
obtains 


T *= 

n 


t+n+1 


g(t) 


P (t) 
R(£+n+l) 


(5.2.21) 


where  P (t)  is  the  (£+n+l)  component  of  P (t) . t cam 
R(£+n+l)  R n 

be  held  constant  by  adjusting  g(t)  at  each  time  instant  to 

account  for  the  time  variations  of  P (t) . This  will  be 

R 

done  in  our  computations  and  will  be  discussed  further  in 
the  next  section. 


f 


I 


Finally  a "motor-noise,"  added  to  the  com- 

memded  control  (see  Fig.  1.1).  Without  this  noise,  the 
model  would  allow  the  operator  to  know  the  precise  value  of 
the  commanded  control  u (t) ; a capediility  the  human  operator 

w 

does  not  possess.  Thus,  the  noise  v_(t)  is  one  source  of 

m 

remnant  which  is  not  referred  back  to  the  model  input. 

Vj||(t)  is  assumed  to  be  zero-meam  white  Gaussian  noise  with 
variance  kernel 

- tj) 

The  magnitude  of  Vjj^(t)  has  been  found  to  vary  with  the  vari- 
ance of  the  control  u(t)  as  follows 

Vm(t)  = p^E{u^(t)}  (5.2.22) 

where  is  a constant  (Refs  5,  23)  . p^^^  has  a typical 

nominal  value  of  .01. 

The  above  discussion  has  summarized  the  mathematical 
description  of  the  optimal  control  model  for  human  opera- 
tors. The  remainder  of  this  chapter  will  be  devoted  to 
discussing  the  solution  to  the  Riccati  equation  and  devel- 
oping equations  for  propagating  mean  states  and  state 
covariances . 


5.3  Solution  to  Riccati  Equation 

In  this  section  we  will  discuss  an  approach  to  solving 
the  steady  state  Riccati  equation  and  some  of  the  charac- 
teristics of  this  equation  as  it  pertains  to  our  problem. 
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Recall  from  Eqs  (5.2.3)  and  (5.2.4)  that  the  first  I states 
of  the  augmented  state  equation  are  uncontrollable.  It 
will  be  shown  below  that  partitioning  the  Riccati  equation 
results  in  the  solution  corresponding  to  the  controllable 
states  being  decoupled  from  the  solution  for  the  uncontrol- 
lable states.  Note  that  it  is  possible  to  formulate  A (t) 
and  a cost  functional  so  that  the  steady  state  Riccati 
equation  will  not  have  a solution.  We  show  in  Appendix  C 
that  our  formulation  results  in  the  existence  of  a solution. 
Now  consider  the  following  partition  of  Eq  (5.2.15)  with 

^R(t)  = 0: 


’ll  ^12  ^lai  r ^ ° 


p T p p 

12  22  23 

p T p T p 

13  23  33 


B + 0 

c 


bT  0 
c 


^11 

^2 

^13 

P T 

P 

P 

12 

22 

23 

P T 

P 

13 

23 

33 

0 0 0 


+ 0 Qd  0 


0 0 0 


11 

12 

13 

T 

P 

p 

12 

22 

23 

T 

P T 

P 

13 

^3 

33. 

1 

g(t) 


I A « 


p 

P 

P 

11 

12 

13 

P T 

P 

P^ 

12 

22 

23 

P T 

P T 

P 

13 

23 

33 

(5.3.1) 
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The  dimensions  of  the  partitions  are  such  that  the  noise 
states  and  system  dyncunic  states  are  separated  (see 
Eq  (5.2.4));  thus 


Ajj  : i X X. 

Aj^  : n X a 

A^  : n X n 

B : n X 1 

c 

Q - : n X n 


Eq  (5.3.1)  leads  to  the  following  set  of  equations: 


^ll^n  ^12^b  ^n^ll  ^b^l2  “ | ^13^13 

* *^12  * *S^22  - ? '’l3''2l  = 0 

*'l2»o  * ^'■X3  * *b’“23  ■ k ''l3''33  ' “ 


= 0 


**22^d  ^6^22  °d  " g ^23^23  “ ° 

P22»c  + *a''23  - ? '’23''33  ' ° 

2P  - — P ^ = 0 

^^23  c g 33  " 


(5.3.2) 

(5.3.3) 

(5.3.4) 

(5.3.5) 

(5.3.6) 

(5.3.7) 


With  g(t)  and  Q,  known,  Eqs  (5.3.5),  (5.3.6),  and  (5.3.7) 
a 

represent  a system  of  [[n(n+l)/2  + n+l]  equations  in  the 
same  number  of  unknowns.  Note  also  that  these  equations 
generate  the  gains  for  the  controllable  states,  independent 
of  the  uncontrollable  states.  Thus  the  optimal  control  for 
the  controllable  states  is  decoupled  from  the  uncontrollable 
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can  be  held  constant  by  the  following  computational 
procedure.  At  t^,  a value  of  is  assumed,  the  Riccati 

equation  is  solved  for  and  then  is  computed  by 

Eq  (5.3.10).  If  the  desired  value  of  is  not  obtained, 
then  q(t^)  is  adjusted  (an  increase  in  g(t)  increases  t^) 
emd  is  recomputed.  This  is  repeated  until  the 

desired  is  obtained  within  the  specified  tolerance.  The 
solution  of  the  Riccati  equation,  P_(t  ) which  gives  the 
desired  value  of  is  used  to  obtain  the  controller  gains 
for  the  interval  t^  to  t^  + At.  At  the  next  computational 
time  (integration  period  At  seconds  later) , ' the  procedure 

is  repeated  with  updated  values  of  A (t) . 

R 

Solving  the  Riccati  equation  requires  a specified 
value  for  the  weighting  matrix  in  the  operator's  cost 
functional.  In  the  tracking  problem  the  operator  is  inter- 
ested in  minimizing  the  mean-squared  tracking  error.  In 
such  a case,  can  be  given  the  form 

■ 0 I o' 

0 I q 

where  q weights  a state  representing  the  difference  between 
the  desired  position  (target)  and  the  actual  position 
(sight).  It  turns  out  that  for  this  class  of  problems, 
holding  constant  results  in  the  controller  gains  being 
Independent  of  the  magnitude  of  q.  This  is  shown  in  Appen- 
dix D.  Thus  specifying  a value  for  results  in  the  ratio 


of  g(t)  to  q being  adjusted  so  that  the  same  control  gains 
are  obtained  for  any  positive  value  of  q.  Thus  we  see  that 
holding  constant  not  only  has  a physical  appeal,  but  is 
also  significant  from  a parameter  identification  standpoint. 
With  this  formulation,  identification  problems  associated 
with  the  controller  gains  are  reduced  to  identifying  the 
proper  value  of  This  completes  the  discussion  of  the 

Riccati  equation  solution;  the  next  section  will  consider 
the  mean  state  and  covariance  propagation. 


5.4  Mecin  State  and  Covariance  Propagation 

We  are  interested  in  using  ensemble  averages  of  simu- 
lator data  to  attempt  to  identify  the  parameters  of  the 
optimal  control  model.  Therefore,  the  mathematical  formu- 
lation must  permit  use  of  the  model  to  predict  the  average 
hximan  operator  response  in  a tracking  problem  in  which 
deterministic  target  motions  are  specified.  That  is  to 
say,  the  target  motion  is  dictated  by  a deterministic  tra- 
jectory that  is  not  known  a priori  by  the  operator  in  the 
tracking  loop.  This  will  be  accomplished  by  using  a tech- 
nique similar  to  Kleinman's  (Refs  20,  23)  and  as  discussed 
by  Jazwinski  (Ref  16) . A determinsitic  component  is 
included  in  the  input  disturbance  w(t)  of  Eq  (5.2.1)  so  that 
the  disturbance  becomes 


w(t)  = w^(t)  + W2(t) 


(5.4.1) 


where  w^(t)  is  zero-mean  white  Gaussian  noise  with 


; -!  ■ s ^ - .-.A  (t  ' V , ' 
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Independent  components  so  that  the  covariance  kernel  is 
diagonal  and  given  by 

E{Wi’'(ti)Wi(t2)}  = Wj^(t^)5(t2  “ t^)  (5.4.2) 

and  W2 (t)  is  a deterministic  disturbance  such  as  target 
acce le  ration . 

From  Eg  (5.4.1)  we  see  that  E{w(t)}  = W2 (t) . As  we 
will  see  below,  this  deterministic  component  will  give  rise 
to  a mean  tracking  error,  which  is  a function  of  time. 

This  mean  response  is  the  result  one  would  expect  to  obtain 
by  ensemble  averaging  the  results  of  many  experimental 
trials  using  the  same  target  trajectory. 

The  covariance  kernel  for  w(t)  is 

E{[w(tj^)  - E{w(t^)>]'^  [w(t^*E{w  (t2>  }]  } 

* W^(t^)6(t^  - t^)  (5.4.3) 

However,  the  operator  has  no  a priori  knowledge  of  w^ (t) . 
Therefore,  he  assumes  (for  lack  of  better  knowledge)  that 
the  disturbance  is  zerd’^mean  white  Gaussian  noise  with  a 
covariance  kernel 

E{w'^(tj^)w(t2)  } = W(tj^)6(t^  - t^)  (5.4.4) 

where  the  diagonal  components  have  values 

W^(ti)  = W^^(t^)  + W2^(t^)  (5.4.5) 
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Thus,  from  the  operator's  viewpoint,  the  deterministic  com- 
ponent contributes  to  the  covariance  of  what  the  operator 
assumes  is  zero-mean  white  Gaussian  noise.  Thus,  the  Kalman 
filter  and  predictor  of  the  model  compute  their  estimates 
as  if  the  disturbance  were  zero-mean  white  Gaussian  noise 
with  a covariance  matrix  W(t)  given  by  Eq  (5.4.4).  This 
approach  results  in  an  adaptive  tracker  (Refs  15,  34)  since, 
as  the  acceleration  of  the  target  is  increased,  the  noise 
disturbance  will  increase  and  thereby  increase  the  gain  of 
the  Kalman  filter.  This  results  in  more  weight  being  given 
to  recent  observations,  which  is  desirable  when  the  target 
is  maneuvering  significantly. 

Kleinman  (Refs  20,  23)  has  modeled  the  plant  disturb- 
ance without  the  white  noise  component;  that  is 

w(t)  = w^  (t) 

where  W2 (t)  is  the  deterministic  disturbance  discussed 
above.  Again,  as  discussed  above,  the  operator  assumes  this 
is  a zero-mean  white  Gaussian  noise.  Kleinman  chose  to  model 
the  variance  kernel  of  what  the  operator  presumes  to  be 
zero-mean  white  Gaussian  noise  by 

w^(t)6(t)  = 2 W2?(t)6(t) 

where  t may  be  viewed  as  a correlation  time  for  w (t) . 
c 2 

The  model  equations  which  are  developed  below  will  assume 
the  white  noise  component  is  present;  however,  we  will  com- 
pare results  obtained  using  the  two  approaches  in  Chapter  VII. 
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As  mentioned  above,  even  though  the  operator  (and 
therefore  the  model  Kalman  filter  and  predictor)  assumes 
the  disturbance  is  zero-mean  white  Gaussian  noise,  we  will 
separate  the  deterministic  and  random  effects  of  the  dis- 
turbance in  developing  the  propagation  equations  for  the 
mean  states  and  the  corresponding  covariances. 

A modified  augmented  state  equation  is  defined  for  the 
purpose  of  deriving  these  equations,  i.e.. 


X (t)  = A (t)x  (t)  + B U (t)  + r w (t) 
a a a a c a a 


(5.4.6) 


where 


Xa(t) 


■x(t)‘ 

Lu(t)_ 


Aa(t)  = 


B = 
a 


r (t)  = 

a 


w (t)  = 

a 


A(t)  I B(t) 

0 I -1/t 


r(t)  I 0 

I 

0 ' 1/T 

I n 

Wj^(t)  + Wj(t) 


v„(t) 

m 


E{w  (t) } = w (t) 

d 4b 


E{[Wa(tj^)  - W2(t^)] 

[Wa(t2>  - W2(t2)]*^} 


»l  1 0 


0 V. 


I 'm 


6(t2  - t^) 


"'-  fi **■> 


Recall  that  U (t)  is  the  commanded  control  given  by 
c 

Eq  (5.2.19). 

Our  objective  is  to  derive  expressions  for  propagating 
the  mean  states 


X (t)  = E{x  (t)}  (5.4.7) 

am  a 

eind  covariance 

P^(t)  = E{[xg(t)  - Xg^(t)]  [Xg(t)  - Xgj^(tf}  (5.4.8) 

Consider  first  the  covariance  of  the  Kalman  filter 
estimate  of  the  delayed  states.  We  will  let  a = t - t where 
T is  the  perceptual  delay  time.  The  Kalman  filter  estimate 

/V 

is  denoted  by  x_ (a ) and  its  covariance  by 

ci 

P (a)  = E{[x  (a)  - X (a)]  [x  (a)  - x (a)]^}  (5.4.9) 

X cl  cL  cl  3 

with  the  augmented  state  vector,  the  operator  observation 
process  is  described  by 

y(o)  = C^(o)x^(a)  + v(t)  (5.4.10) 

where 

C^(a)  = [C(t)  I D(t)]  (5.4.11) 

and  v(t)  is  zero-mean  white  Gaussian  noise  with  independent 
components  and  variance  kernels 

E{v^(t^)v^(t2)>  = V^(tj^)6(t^  - t^) 
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Then  the  covariance  of  the  Kalman  filter  estimate  is  propa- 
gated by 

P (a)  = A (o)P  (a)  + P (a)A'^(a)  + r (a)W(t) (a) 

1 a 1 la  a a 

- P (a)c'‘^(a)v“^(cr)C  (a)P  (a)  (5.4.12) 

X d d X 

Also,  the  estimate  of  the  delayed  states  obeys 
X 

x(<J)  = A^(a)x  (a)  + B U (a) 
a a a c 

+ P (a)C^(a)v"^(a)  [y(a)  - C (a)x^(a)] 

X a Si  cL 

(5.4.13) 

Subtracting  Eq  (5.4.13)  from  Eq  (5.4.6),  we  obtain  the 

equation  for  propagating  the  estimation  error  of  the  delayed 
A 

states,  e- (a)  = x (a)  - x (a) 

1 a a 

e (a)  = A (a)e,  (a)  - P,  (a)c'’’^(0)v“^  (a)[c  (a)e  fa)  + v(a)] 
j.  a j.  j.  a a 1 

+ r (a)W(a)  (5.4.14) 

a 

Making  the  definitions 

G(0)  = P,  (0)C^(0)V"^(0) 

X d 

A^(0)  = A (0)  - G(0)C  (0)  (5.4.15) 

fa  a 

where  G(0)  is  the  estimator  gain  and  A^(o)  is  the  closed- 
loop  filter  gain.  Eq  (5.4.14)  becomes 


e,  (o)  = A (o)e.  (0)  - G(0)v(0)  + r (0)W(0) 


(5.4.16) 


If  we  further  define  the  closed  loop  control  matrix  by 


A (a)  = A (o)  - B X (a) 
c a a c 


where  X . is  given  by  Eq  (5.2.19)  for  i=l,  2,  •••  i+n 
2Uid  X^  £,+n+l  = 0*  expression  for  propagating  the 

estimate  of  the  delayed  states  (a)  in  Eq  (5.4.13)  becomes 

Si 

• 

X (a)  = [a  (a)  - B X (a)]x  (a)  + G(a)[c  (a)e  (a)  + v(a)] 
a a aca  al 

A 

= A (a)x^(a)  + G(a)[c  (a)e  (a)  + v(a)]  (5.4.17) 

C a cl  X 


Now  the  predictor  generates  the  best  estimate  x (t) 
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A 

from  the  Kalman  filter  output  x (a).  This  is  obtained  by 

Si 

evaluating 


X (tjo)  = E{x  (t)  I X (a)} 


(5.4.18) 


This  leads  to  the  following  equation  for  propagating 

A 

X (t  |a)  ; 

Ci 

x^(t|a)  = A^(t)x^(t|a)  + B^U^(t) 

**<|)(t,  t-T)G(a)[c  (a)e  (a)  + v(a)] 

cl  X 

(5.4.19) 

Now  the  transition  matrix  <|)(t,  t-r)  obeys 


3<t>(ty  t-T) 

3t 


= A^(t)4»(t,  t-T) 
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For  systems  where  A (t)  changes  only  slightly  between  t-x 

a 

and  t (recall  x « .2  sec),  then 

♦ (t,  t-x)  3 e 

Also  making  the  definition 

K(a)  = e G(a) 

Substituting  these  quantities  into  Eq  (5.4.19)  yields 
x^(t,|a)  = A^(t)x^(t|a)  + K(a)[c^(a)ej^(a)  + v(a)]  (5.4.20) 

The  prediction  error,  error  introduced  in 

A I A 

estimating  x_(t|0)  given  x_ (a)  and  is  obtained  by 

d ci 

t Aa(t-e) 

e (t)  = / e r,(e)w(e)de  (5.4.21) 

^ t-x  ® 

Combining  the  above  results  yields  the  total  state, 

Xjj  (t) , as 

Xa(t)  = x^(t/a)  + e^®^‘^^^ej^(a)  + C2(t)  (5.4.22) 

The  mean  state,  x,„(t),  is  determined  by  evaluating 

cuu 

the  expected  value  of  the  three  terms  of  Eq  (5.4.22);  i.e., 

.s  A^(a)x 

= E{x_(t))  + e E(e,  (a)}  + E{e,(t)} 

am  a x z 

(5.4.23) 

From  Eqs  (5.4.20),  (5.4.16),  and  (5.4.21)  we  can  determine 
the  following  equations  for  propagating  the  expected  values. 
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AlikyiUk 


X (tia)  = A (t)x  (t|a)  + K(a)C  (a)e  (a)  (5.4.24) 

&fii  I c Am  ^ J fn 


[-r] 


. Im 


Wj(o) 


(5.4.2S) 


‘ rr(e)l 


(5.4.26) 


where 


X (t,la)  = E{x  (tja)> 
am  a 


e (a)  = E{ie  (a)} 
Im  1 


le,  (t)  = E{e  (t)} 
2m  2 


Note  from  the  above  coupled  set  of  matrix  equations  that 
the  me  cm  states  occur  because  of  the  deterministic  quantity, 
w (t) , in  the  disturbance. 

A 

Now  the  following  definitions  are  established  for  cal- 
culating covariance  and  crosscovariance  terms. 


“ E{[e^(a)  - 

(c)  = E{[x  (t)  - X (t)][e,(o)  - e (a)]*^} 
J a am  X xm 

^ (t)  = E{[x  (t)  - x^„(t)][x  (t)  - x„(t)]*^} 
4 a am  a am 


(5.4.27) 


. 

J, -V  • X vr^,. 


From  these  definitions  and  the  previous  results,  the 
following  equations  are  obtained  for  propagating  the  above 


quantities : 


® A^(a)P2(a)  + P2(o)A^(o)  + ^^(a) 


W^(a)  I 0 


0 1 V„(cr) 

* lU 


T 

r (a) 
a 


+ G(a)V(a)G  (a) 


(5.4.28) 


PjCt) 


P4(t) 


A^(t)P2(t)  + P2(t)A^(o)  + K(a)C^(a)[P2(o)  - Pj^(a)] 

(5.4.29) 

A (t)P,  (t)  + Pjt)A'^(t)  + K(a)C  (o)p'J(t) 
c 4 4 c a .3 


+ P,  (t)c'^(o)K'^(a)  + K(a)V(o)K‘^(o) 
3 a 


(5.4.30) 


Pg(t) 


t A.(t-e) 

/ e r (e) 

t-T  ^ 


W,  (e) » 0 1 


0 'V„(e) 

I m 


A^(t-e) 

r^(e)e  de 

cl 


(5.4.31) 


Finally  the  covariance  of  the  total  states,  P^(t),  is 
determined  by 


p (t)  = e 

d 


A_(o)t  A^^)t  A^(o)t 

P2<^)e 


A»(a)T 


+ e 


pj(t)  + P4(t)  + Pg(t) 


(5.4.32) 
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The  mean  and  covariance  of  the  operator-observed  data 

are 

E{y(t)}  = y^(t)  = C (t)x_(t)  (5.4.33) 

m a am 

E{[y(t)  - y (t)][y(t)  - y^Ct)]"^}  = C (t)P  {t)c’’(t) 
m m a a a 

(5.4.34) 

In  this  chapter  we  have  described  the  optimal  control 
model  for  human  response  in  some  detail.  The  model  has 
several  parameters,  the  values  of  which  influence  the 
response  of  the  system  in  simulating  the  human  operator. 
These  parameters  are 

Tjj  = neuromotor  delay  (Eq  (5.2.16)) 

T = perception  time  delay  (Eq  (5.2.6)) 

= proportionality  constants  for  calculating  obser- 
vation noise  variance  kernels  (Eq  (5.2.8)) 

0 = proportionality  constant  for  calculating  motor 

m 

noise  variance  kernel  (Eq  (5.2.22)) 

= variance  of  white  noise  component  of  plant 
disturbance  (Eq  (5.4.2)). 

These  parameters  appear  in  the  coupled  sets  of  matrix 
differential  equations  which  propagate  the  mean  states  and 
covariances  of  the  model  response.  This  response  would 
correspond  to  ensemble  averaging  many  trials  of  a human 
operator  in  a tracking  situation.  The  development  has  the 


i 


advantage  that  one  can  predict  the  mean  response  to 
deterministic  target  trajectories.  In  the  next  chapter 
we  will  describe  a simulator  system  which  was  used  to 
obtain  data  of  operator  response  in  tracking  maneuvering 
targets  flying  deterministic  trajectories.  This  will  set 
the  stage  for  a nontrivial  parameter  identification 
problem  as  discussed  in  Chapters  II,  III,  and  IV. 


Chapter  VI 


SIMULATOR  SYSTEM  DYNAMICS,  STATE  SPACE 
MODELING,  AND  RESULTS  OF  EXPERIMENTAL  TRIALS 

The  purpose  of  this  chapter  is  to  describe  the  simu- 
lator system  which  was  used  to  obtain  experimental  data 
with  operators  performing  a tracking  task.  The  simulator 
dynamics  will  be  described  and  the  state  space  equations 
will  be  developed  to  correspond  with  the  development  of 
the  previous  chapter.  In  addition,  visual  thresholds  will 
be  discussed  in  the  context  of  the  simulator  display,  as 
well  as  an  indifference  threshold  effect  due  to  the  finite 
size  of  the  displayed  target.  The  deterministic  target  tra- 
jectories which  were  used  for  the  simulator  trials  and  the 
data  resulting  from  ensemble  averaging  many  trial  runs  will 
be  presented  in  this  chapter  as  well. 

6.1  General  Description  of  Simulator 

The  simulator  used  in  this  research  is  maintained  and 
operated  by  the  Aerospace  Medical  Research  Laboratory, 
Wright-Patterson  AFB,  Ohio.  Only  the  portions  of  the  system 
pertinent  to  the  discussion  will  be  presented. 

The  system  uses  two  people  in  the  tracking  scheme. 

Each  person  has  an  independent  set  of  controls  for  tracking 
one  axis  (azimuth  or  elevation) . The  operators  receive 
tracking  information  by  way  of  a television  monitor.  A 
large  screen  presents  data  only  for  the  first  ten  seconds 
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to  allow  acquisition.  After  that  both  operators  receive 
data  from  a small  4.5-inch  by  4.5-inch  monitor.  The 
sighting  reticle  consists  of  pairs  of  narrow  vertical  and 
horizontal  lines  with  a spacing  of  about  .002  inches. 

Each  operator  is  provided  with  handwheels  to  position  the 
geometrical  center  of  the  TV  camera  which  provides  the 
visual  signal  to  the  display  monitor.  The  elevation  oper- 
ator is  attempting  to  maintain  the  horizontal  lines  on  the 
target  and  the  azimuth  operator  controls  the  vertical  pair 
of  lines.  The  operators  have  the  choice  of  either  position 
control  or  a rate-aided  mode  of  control  (see  Section  6.2 
on  simulator  dynamics) . Only  one  mode  select  switch  is 
provided  so  that  each  axis  is  always  tracked  in  the  scune 
mode.  The  usual  operation  is  to  acquire  the  target  using 
the  position  control  followed  by  operation  in  the  rate- 
aided  mode  while  tracking. 

6.2  Simulator  Dynamics 

The  geometry  of  the  tracking  task  is  shown  in  Fig.  6.1. 
From  this  figure  we  see  that 

0^  = azimuth  angle  of  target 
0g  = elevation  angle  of  target 
0g^  = azimuth  angle  of  sight/camera 
0gg  = elevation  angle  of  sight/camera 


% 
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Fig.  6.1  Tracking  Geometry 


The  operators  view  the  traverse  error  (in  the  slant  plane) 
and  the  elevation  error  on  the  monitor;  therefore  we  make 
the  following  additional  definitions : 

Gq,  = traverse  angle  of  target  in  slant  plane 
0g^  = traverse  angle  of  sight /camera 

so  that 

®r  = Sa  =°=®e  (6.2.1) 

®ST  = ®SA 


Prom  Fig.  6.2,  the  following  transfer  functions  are 
obtained: 

For  position  control 


0 

U 


SA 

PA 


2 

s 


+ 6s  + a 


(6.2.3) 


0 


SE 


K, 


U 


PE 


s + 6s  + a 


(6.2.4) 


For  rate  control 


and 


0 

SA 


K 


s(s  + 6s  + a) 


0 

SE 


K 

4 


+ 6s  + a) 


(6.2.5) 


(6.2.6) 


The  following  values  for  the  transfer  function  constants 
correspond  to  those  of  the  simulator  used  for  obtaining 
experimental  data: 


o = 64 

II 

r— 1 

3.2 

Kj  = 3.156 

6 = 12 

K3  = 

.8 

K4  = 8.533 

6 • 3 Simulator  State  Space  Equations 

The  state  space  equations  will  be  developed  to  corre- 
spond with  the  optimal  control  model  presented  in  Chapter  V. 
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The  states  of  the  disturbance  noise  and  system  dynamics  obey 
the  equation 


x(t)  = Ax(t)  + r w(t)  + B(t)u(t) 


(6.3.1) 


where,  assuming  rate-aided  control. 


x(t)  = 


0 - 0 
ST  T 


0 - 0 
SE  E 


0 0 0 0 

0 0 10 

0 -a  -B  0 

-l.  'l  0 0 


0 0 C 

0 0 ( 

0 0 C 

0 0 -] 


0 0 
0 0 
0 0 


1 0 


0 -a  -B  0 


w(t)  = (t)  + (t) 


B(t)  = 


K2  cosQg(t) 


uTt)  * 


As  discussed  in  Section  5.4,  Wj^(t)  is  zero-mean  white 
Gaussian  noise  and  ”2^^^  ® deterministic  disturbance. 

We  will  use  target  angular  acceleration  as  the  deterministic 
portion  of  the  disturbance  so  that 


similarly r the  elevation  axis  can  be  described  by  an 


equation  of  the  form  of  Eq  (€.3.1),  with 


x(t)  = 


e - e 

SE  E 


0 0 


1 0 


-e  0 


U(t)  = Uj^(t) 


0 0 


w(t)  = w (t)  + 0 (t) 
1 E 


Note  that  B is  not  a function  of  time  here.  The  state  vec- 


tor in  each  of  these  decoupled  equations  has  the  first  ele- 


ment as  state  for  the  noise  model  and  the  remaining  three 


are  the  states  of  the  system  dyneimics  (see  Eq  (5.2.3)}. 


The  variables  for  each  axis  displayed  to  the  operators 


yd<t)  = C x(t) 


(6.3.2) 


where 


[ 0 0 0 l“ 

-110  0. 


That  is,  each  operator  is  assumed  to  observe  angular  error 


and  angular  error  rate.  As  discussed  in  Chapter  V,  this  is 
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consistent  with  the  usual  assiimption  that  if  a variable  is 
displayed  visually  to  the  operator,  then  the  rate  of  change 
of  that  quantity  is  also  discerned. 

Each  operator  wemts  to  miniitiize  his  own  cost  functional, 
i.e.,  for  azimuth 

+ g(t)Uj^(t)]dt} 

0 ' 

0 
0 

q . 

The  elevation  operator  has  a similar  cost  functional. 

6,4  Simulator  Riccati  Equation  (see  Section  5.3) 

The  remaining  discussion  will  be  confined  to  the  azi- 
muth axis,  with  obvious  extension  to  the  elevation  axis. 

We  will  approach  the  solution  to  the  Riccati  equation  for 
the  simulator  as  discussed  in  Section  5.3.  An  augmented 
state  vector  is  defined  as  in  Eq  (5.2.11) 


= lim  E {J-  / [x^(t)Q  x(t) 

t--*^  t. 


where 


Q = 


0 0 0 
0 0 0 
0 0 0 
0 0 0 


XR(t)  = 


x(t) 

u(t) 


with 


where 


XR(t)  - Aj^^t)Xj^(t)  + Bj^u(t) 


A^(t)  = 


A B(t) 

0 0 


(6.4.1) 


(6.4.2) 


and 


L 


Note  that  the  subscript  R is  used  to  denote  the  augmented 
state  equation  leading  to  the  Riccati  equation  as  distin- 
guished from  another  augmented  state  equation  in  Section  6.5. 
The  optimal  gains,  X,  are  obtained  by 


X(t)  =1  Pj^(t)  (6.4.3) 

Recall  that  P_(t)  is  the  solution  of  the  Riccati  equation 

Xx 

at  time  t with  t^^“  ; i.e.,  the  quasi-static  steady  state 
solution  to 


Pc<t)B„ 


g(t) 


+ Q =0 
R 

(6.4.4) 


where 


Partitioning  the  Riccati  equation  as  shown  in  Eq  (5.3.1) 
leads  to  the  following  definitions  for  the  simulator  case: 


0 

0 

-1 
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1 

0 

iH 

0 

1  

0 

II 

1 

0 0 

ca  0 

1 

0 pH 

1 

1 

B (t)  = 
c 

K,  cos0  (t) 
2 E 

0 

Note  that  time  variations  occur  in  the  azimuth  case  due  to 
the  crosscoupling  of  the  elevation  angle  (appearing  in 

that  do  not  occur  in  the  elevation  axis.  Then  the 
equations  corresponding  to  Eqs  (5.3.2)  through  (5.3.7) 
become 


12 


■>•[0  0 


(6.4.5) 


12 


0 10 
-o  -6  0 

10  0 


+ [o  0 -l]Pjj  - I - 0 (6.4.6) 


12 


Kj  COS0j.(t) 
0 


*[0  0 -1>33  - 5 W33  = ® 


' 0 1 o' 

1 

0 

1 
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M 

1 

'0 

0 

0 

1  

®*22 

-a  -8  0 

+ 

1-8  0 

**22  * 

0 0 0 

M 

0 

0 

— 1 

0 

0 

0 

i 

_0  0 q 

1 T 

-ip  P = 0 
g 23  23 


(6.4.8) 


0 

0 

-0 

l‘ 

^22 

Kj  cos©  (t) 

+ 

1 

0 

0 

0 

0 

0 

23 


ip  p 

g 23  33 


* 0 


(6.4.9) 
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0 


*23  "^2  ’ g '*33  = ° 

0 


(6.4.10) 


Also  recall  from  Section  5.3  that  if  = g(t)/P22(t)  is 
held  constant,  this  completely  determines  the  gains,  X (t) , 
regardless  of  the  value  of  q.  The  procedure  used  here  was 
to  specify  a value  of  set  q = 1,  and  then  adjust  g(t) 
to  obtain  the  desired  value  of  =;  g(t)/P22(t).  (A 'value 
within  ± .001  of  the  desired  value  was  considered  adequate.) 
With  q and  g(t)  specified,  Eqs  (6.4.5)  through  (6.4.10)  can 
be  solved  for  the  required  unknowns  to  obtain  the  feedback 
gains  X (t) . The  details  of  a method  for  obtaining  a 
solution  are  discussed  in  Appendix  E.  As  shown  in 
Eq  (5.3.9),  X (t)  is  determined  by 


=iW[^3  : -23  1 -33] 

6.5 

The  augmented  state  Eq  (5.4.6)  takes  the  following 
form  for  the  simulator  azimuth  axis: 

x^(t)  = A (t)x  (t)  + B U (t)  + r w (t)  (6.5.1) 

d 3 d 3 Q 3 3 


simulator  Mean  State  and  Covariance  Propagation 
(see  Section  5.4) 


’x(t) 

,u(t) 


0 0 0 0 

0 10  0 

a -B  0 K2  cos0g(t) 

10  0 0 

0 0 0 -lAn 


0«(t) 


X^(t)x(t) 


with 


Ki  = T X,  (t)  i = 1,  2,  3,  4 

cl  n 1 

With  the  sJaove,  the  following  equations  result  for  propa- 
gating the  mean  states  as  derived  in  Section  5.4: 

"‘am**’  ' ® 

where  x (t)  = E{x  (t) } and  ct  = t - x. 

cUKl  cl 

A 

The  equations  for  propagating 

are  expressed  in  Eqs  (5.4.24)  through  (5.4.26). 

The  covariance 

Pa(t)  = E{[x^  - x^(t)]  [x^(t)  - Xg^(t)]  (6.5.3) 


is  determined  from 


P^(t)  = e 


Aj^(a)x 


P2(tJ)e 


A^(a)x 


+ P3(t)  e 


A^(a)x 

a 


A- (a) XT 


P3(t)  + P^(t)  + P5(t) 


(6.5.4) 


where  P^(t),  P^it) , P^(t) , P^(t),  and  P^ (t)  are  propagated 
by  Eqs  (5.4.12)  and  (5.4.28)  through  (5.4.31). 


6.6  Threshold  Effects 

As  previously  discussed,  the  simulator  operators 
observe 


y(t)  = . 


®ST  " ®T 


©ST  “ ®T 


' . ■ ■ 5:,^ 
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I, 
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nominal  position.  A survey  was  taken  of  seven  tracking 
teams  (14  people)^  and  it  was  foxind  that  the  average  region 
of  indifference  was  about  ±8  feet  in  the  traverse  plane  and 
±5  feet  in  the  elevation  plane.  Thus  a second  threshold  is 
introduced  which  is  a function  of  target  range;  i.e.,  the 
magnitude  of  the  threshold  increases  as  the  size  of  the 
image  on  the  screen  increases.  Incorporation  of  this  effect 
is  considered  a significant  refinement  since  this  same 
apparent  phenomena  has  been  observed  in  other  experimental 
data,  but  not  adequately  explained.  In  view  of  this,  a 
threshold  was  Incorporated  as  follows: 


Position  Threshold  = 


••3  Ax 

Max  1.175  X 10  rad. 


[• 


Range 


Rate  Threshold  = 


-3 


Max  .628  x 10  rad,  .628  x 10 


_3  r Pos.  Thres, 


■P 

_ • 


175  X 10 


-3 


where 

Ax  = 8 feet  for  tr^erse  axis 
Ax  = 5 feet  for  elevation  axis 


^Eight  teams  of  two  people  each  were  used  in  the  simu- 
lation runs;  one  team  was  not  available  for  this  survey. 
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The  threshold  effect  is  incorporated  into  the  model  as 
in  Ref  2.  Each  perceived  output  is  modified  according  to 


where  the  threshold  nonlinear  function  is 


(6.6.1) 


f (X)  = 


X - 


X + a. 
1 


X > a. 
— 1 


“ < X < a^ 


X < - a. 
— 1 


(6.6.2) 


To  facilitate  accounting  for  this  nonlinear  effect  in  the 

A 

propagation  equations , we  want  to  find  a function  f (y)  such 
that  the  difference 


d(t)  = f(y(t))  - f(y)  . y(t) 


(6.6.3) 


is  minimized  in  the  mean  squared  statistical  sense.  When 
y(t)  is  assumed  to  be  a zero-mean  Gaussian  random  variable 
with  variance  o,  it  has  been  shown  (Ref  21)  that 


f(y)  = erfc 


a 


\f2 


(6.6.4) 


where 


erfc(b)  = 1.  - 


, b 2 
_2_  -w 

j—  ; p dw. 
1 TT  o 


Thus  it  is  assumed  that  the  operator  perceives 

y(t)  = f(y)y(t)  + v(t) 


(6.6.5) 
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Recall  from  Eq  (5.4.10)  that 


y(t)  = C X (t)  + v(t) 

a cl 

where 

so  that,  considering  threshold  effects 

y(t)  = f(y)C  X (t)  + v(t) 

d Si 

A 

Now  define  C'  = f(y)C  , to  obtain 

cl 


(6.6.6) 


(6.6.7) 


y(t)  = C X (t)  + v(t)  (6.6.8) 

a 

Recall  that  v(t)  is  zero-mean  white  Gaussian  noise  with 
independent  components  and  variance  kernels  (t)  . Since 

fp 

only  the  quantity  C_  V (t)  C appears  in  the  model 

oL  cl 

equations,  an  equivalent  result  is  obtained  by  assuming 
that  the  perceived  data  is 


t-1 


y(t)  = y^(t)  + v(t)  f (y) 


= C X (t)  + v(t)  f“^(y) 
a a 


(6.6.9) 


The  observation  noise  covariance  without  threshold  effects 
is  obtained  by 


var(y) 


(6.6.10) 
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Considering  threshold  effects,  the  covariance  is  obtained 
by 


Vi  = Cp^  var(y)] 


''-2 

f (y) 


(6.6.11) 


where 


f(y) 


erfc 


6.7  Target  Trajectories 

Two  target  trajectories  are  used  to  generate  simulator 
data.  Trajectory  1 is  a fly-by  (see  Fig.  6.7.1).  The  tra- 
jectory initiates  at  = 10000  feet,  y^  = 5000  feet,  and 
Zq  = 5000  feet  and  flies  straight  and  level  as  shown  at  a 
velocity  of  500  feet/sec.  Trajectory  2 is  shown  in 
Fig.  6.7.2  with  coordinates  given  in  Table  VI -I.  Equations 
used  to  calculate  angular  position,  velocity,  and  accel- 
eration for  the  fly-by  trajectory  are: 


e. 


0 

T 


• • 

e_ 


-1  r^o  ~ ^o^~ 

L 

COS0 

£ 

• * 

0 COS0  - 0 0 sin0 

A E A E E 

0 COS0  - 0-  0„  sin0„  - 0-  sin0 
A E^E  EAE  E 

" ®A  ®*E  " ®A  ®E^ 


0. 


COS0„  = 
E 


tan 


-1 


- ^o^ 


(6.7.1) 

(6.7.2) 

(6.7.3) 


(6.7.4) 
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Table  6.1 


TIME  VS.  POSITION  (TRAJECTORY  2) 


Time 

- X 

y 

z 

.00 

-11318 

8723 

7441 

.50 

-11129 

8389 

7441 

1.00 

-10900 

8455 

7441 

1.50 

-10691 

8321 

7441 

2.00 

-10482 

8187 

7441 

2.50 

-10273 

8053 

7441 

3.00 

-10064 

7919 

7441 

3.50 

-9855 

7785 

7441 

4.00 

-9645 

7651 

7441 

4.50 

-9437 

7517 

7441 

5.00 

-9228 

7383 

7441 

5.50 

-9219 

7249 

7441 

6.00 

-8810 

7115 

7441 

6.50 

-8601 

6981 

7441 

7.00 

-8392 

6847 

7441 

7.50 

-8183 

6713 

7441 

8.00 

-7974 

6579 

7441 

8.50 

-7763 

6445 

7441 

9.00 

-7538 

6311 

7441 

9.50 

-7347 

6177 

7441 

10.00 

-7133 

6043 

7441 

10.50 

-6929 

5909 

7427 

11.00 

-6720 

5775 

7378 

11.50 

-6511 

5641 

7285 

12.00 

-6302 

5507 

7149 

12.50 

-6093 

5373 

6969 

13.00 

-5884 

5239 

6745 

13.50 

-5675 

5103 

6544 

14.00 

-5452 

4098 

6336 

14.50 

-5219 

4911 

6122 

15.00 

-4980 

4839 

5902 

15.50 

-4736 

4778 

5676 

16.00 

-4487 

4726 

5445 

16.50 

-4235 

4879 

5210 

17.00 

-3978 

4837 

4972 

17.50 

-3719 

4598 

4729 

18.00 

-3455 

4561 

4483 

18.50 

-3189 

4526 

4234 

19.00 

-2919 

4490 

3981 

19.50 

-2644 

4458 

3727 

20.00 

-2361 

4423 

3478 

20.50 

-2065 

4388 

3236 

21.00 

-1752 

4350 

3013 

21.50 

-1425 

4311 

2811 

22.00 

-1084 

4270 

2630 

22.50 

-732 

4228 

2471 

Table  6.1  (continued) 


Time 

- X 

y 

z 

23.00 

-370 

4184 

2336 

23.50 

-1 

4140 

2223 

24.00 

374 

4095 

2135 

24.50 

752 

4050 

2070 

25.00 

1132 

4004 

2030 

25.50 

1511 

3059 

2014 

26.00 

1888 

3014 

2022 

26.50 

2260 

3869 

2053 

27.00 

2627 

3829 

2110 

27.50 

2990 

3812 

2171 

28.00 

3331 

3826 

2242 

28.50 

3706 

3873 

2302 

29.00 

4064 

3951 

2554 

29.50 

4392 

4039 

2400 

30.00 

4713 

4198 

2441 

30.50 

5024 

4364 

2477 

31.00 

5313 

4530 

2509 

31.50 

5802 

4752 

2541 

32.00 

5891 

4946 

2573 

32.50 

6180 

5140 

2625 

33.00 

6469 

5334 

2637 

33.50 

6758 

5328 

2659 

34.00 

7047 

5722 

2701 

34.50 

7336 

5816 

2733 

35.00 

7625 

6110 

2765 

35.50 

7914 

6304 

2797 

36.00 

8203 

6498 

2829 

36.50 

8492 

6692 

2861 

37.00 

8781 

6886 

2893 

37.50 

9070 

7280 

2926 

38.00 

9359 

7274 

2957 

38.50 

9648 

7468 

2939 

39.00 

993^ 

7662 

3021 

39.50 

10226 

7856 

3053 

40.00 

10515 

8050 

3085 

40.50 

10804 

8244 

3117 

41.00 

11093 

8438 

3149 

41.50 

11382 

8532 

3181 

42.00 

11671 

8826 

3213 

42.50 

11960 

9020 

3245 

43.00 

12249 

9214 

3277 

43.50 

12538 

9408 

3309 

44.00 

12827 

9602 

3341 

44.50 

13116 

9796 

3373 

45.00 

13405 

9990 

3405 
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6.8  Simulator  Data 

Eight  teams  of  two  people  each  were  used  in  the  simu- 
lation runs.  The  tracking  error  data  was  recorded  as  a 
function  of  time  on  magnetic  tape  for  many  runs  of  each 
te2un  for  both  of  the  trajectories.  Ensemble  averages  were 
computed  at  .2  second  intervals  over  56  runs  for  each  tra- 
jectory. In  other  words,  at  each  time  point,  the  mean 
tracking  error  was  computed  for  each  axis  as  follows 


1 N 

®sm(tj>  = k I 


i=l 


(6.8.1) 


where 


0 (t . ) = mean  tracking  error  at  time  t . 
sm  j ’ ] 

0g^  (t j ) = tracking  error  at  time  t ^ for  run  i 
N = total  number  of  runs  = 56 


The  standard  deviation  was  computed  using  the  relation 


fl  N 2 2 1 2 


(6.8.2) 


where 


o (t. ) = standard  deviation  at  time  t. 
s j D 


No  attempt  was  made  to  subgroup  data  by  teams  since  the 
objective  here  was  to  obtain  a large  number  of  trials  over 
a crosssection  of  teams. 
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Plots  of  this  data  for  each  axis  for  the  two  trajec- 
tories are  shown  in  Figs.  6.8.1,  6.8.2,  6.8.3,  and  6.8.4. 

The  large  errors  initially  occur  during  acquisition  of  the 
target  by  the  trackers.  On  each  of  the  figures  we  have 
denoted  some  important  time  points  associated  with  the  tra- 
jectories. For  example,  for  the  fly-by  of  Trajectory  1, 
the  time  of  closest  approach  (crossover)  is  indicated.  Note 
that  the  point  of  crossover  is  not  so  close  that  track  is 
lost.  For  Trajectory  2,  the  times  of  azimuth  crossover  and 
minimum  altitude  are  indicated. 

In  the  next  chapter  we  will  use  this  data  in  con- 
junction with  a gradient  computational  procedure  to  eval- 
uate the  optimal  control  model  parameters. 


no.  6.8.1 

AZIMUTH  TRACKING  ERROR 
TRAJECTORY  ONE 


TRACKING  ERROR 

(e.  - 


MEAN  ERROR 


CROSSOVER  (X=0) 


TIME  (SEC) 


30.00 


Chapter  VII 


IDENTIFIABILITY  OF  OPTIMAL  CONTROL  MODEL 
AND  COMPUTATIONAL  RESULTS 


In  this  chapter,  the  identifiability  of  the  optimal 
control  model  parameters  will  be  evaluated,  and  the  results 
of  computations  to  estimate  the  values  of  the  parameters 
will  be  presented.  With  the  structure  of  the  model  for 
hiiman  performance  assumed  to  be  that  presented  in  Chapters 
V and  VI,  there  are  several  parameters  whose  values  need  to 
be  estimated  based  on  measured  data.  The  parameters  of 
interest  are; 

T neuromotor  delay  time 

T -—perception  time  delay 

pj^*—proportionality  constants  for  calculating  obser- 
vation noise  variance  kernels 

Ppi — proportionality  constant  for  calculating  motor 
noise  variance  kernel 

Wi-  -variance  of  white  noise  component  of  scalar  plant 
disturbance. 

These  parameters  are  viSWed  as  elements  of  a vector  (p;  i.e.. 


<!>  = [tn  T p^  P2  W^]'^ 


(7.1) 


Recall  from  Chapter  V that 

Xam(t)  = mean  states  of  optimal  control  model 

Pg (t)  = covariances  of  states  of  optimal  control 
model 


As  discussed  previously,  one  can  construct  a set  of  meas- 
ured data  by  computing  ensemble  averages  and  variances  of 
many  trials  on  a simulator.  An  example  of  this  type  of 
data  is  shovm  in  Figs.  6.8.1  through  6.8.4.  Suppose  we 
define  a vector  P (t) , which  is  formed  from  the  diagonal 

ciXX 

components  of  P (t) . With  this  definition  we  can  form  a 

Ck 

state  vector 


x(t) 


eon 


(t) 


axx 


(t) 


(7.2) 


Then  the  coupled  set  of  vector  matrix  differential 
equations  for  propagating  the  mean  states,  x,„(t),  and  the 

alu 

covariances  P^(t)  derived  in  Chapter  V can  be  used  to  form 
the  nonlinear  differential  equation 


x(t)  = f(t,  x(t),  W2(t),  (|))  (7.3) 

where  x(t)  and  (|)  are  as  defined  above  and  W2  (t)  is  the 
deterministic  component  of  the  input  disturbance  defined  in 
Eq  (5.4.1).  In  the  context  of  Eq  (7.2),  W2 (t)  is  regarded 
as  the  input  or  control  signal.  Note  that  the  control  out- 
put of  the  human  model,  u(t),  is  one  of  the  components  of 
the  state  vector  x^^(t)  and  thereby  its  average  value  and 
standard  deviation  are  components  of  x(t)  in  Eq  (7.3).  The 
observations  are  viewed  as  a function  of  the  states  x(t); 
i « e . , 


y(tj^)  = h(t^,  x(t) , 4>) 


f 


•jr»*  ~ 

y.  V 


From  the  simulator  operation,  one  obtains  data,  z(t^), 
vhich  are  the  measured  counterparts  to  the  model  values 
y(t^)  . The  above  discussion  is  intended  to  put  the  optimal 
control  model  parameter  identification  problem  into  the 
context  of  identifiability  of  a nonlinear  system  as  dis- 
cussed in  Chapters  II  and  III.  In  the  next  two  sections  we 
will  use  the  results  of  these  chapters  to  examine  the  local 
and  global  identifiabi J 5 ty  of  the  optimal  control  model 
parameters. 


7.1  Local  Identifiability  of  the  Optimal  Control  Model 
Recall  from  Theorem  3.2  that  the  parameter  vector  (p 
is  locally  identifiable  under  the  observation  process  if 


k 

M((|))  = I 

i=l 


■9y(ti)l^  1 

94)  J 

L 

(7.1.1) 


is  nonsingular. 

The  matrix  M(({>)  will  be  computed  in  conjunction  with  the 
gradient  computations  used  to  estimate  the  model  parameter 
values.  In  Section  7.4  we  will  see  that,  for  each  of  the 
cases  where  computatiortS’  were  made,  the  determinant  of  M(<()) 
was  nonzero.  If  the  determinant  of  M((l))  is  nonzero,  then 
M((}))  is  nonsingular,  thereby  indicating  that  the  parameter 
vector  <j)  of  the  optimal  control  model  is  locally  identifi- 
able under  the  observation  process  used  for  the  computations. 
Recall  also  from  Theorem  4.8  that  if  we  find  a critical 
point,  4i*,  of  the  least  squares  cost  functional  J((|))  and 
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is  nonsingular,  then  this  not  only  assures  local 
Identifiability,  but  also  assures  that  (|)*  is  a local  mini- 
mizer  of  J(<ti).  This  will  be  discussed  further  in  a later 
section  when  computational  results  are  presented. 

7.2  Global  Identif iability  of  the  Optimal  Control  Model 
In  Chapters  II  and  III,  sufficient  conditions  were 
developed  to  establish  the  global  identifiability  of  the 
parameter  vector  (}>.  The  complexity  of  the  equations  of  the 
optimal  control  model  makes  it  difficult  to  use  these 
theorems;  however,  we  will  use  some  of  the  results  of 
Chapter  III  to  draw  certain  conclusions  concerning  the 
global  identifiability  of  4>.  In  the  following  sections  we 
will  examine  the  global  identifiability  of  each  of  the 
pareuneters  of  the  optimal  control  model. 

7.2,1  Global  Identifiability  of 

Consider  an  observation  of  the  form 

y(t)  = H x^jj^(t)  (7.2.1) 

where  x (t)  is  the  vector  of  mean  states  given  by 
am 

Eg  (5.4.22),  i.e., 

^ . A-(o)t 

X (t)  = X (t  o)  + e e (o)  + e (t)  (7.2.2) 

am  am  im  2m 

We  will  apply  Theorem  3.9  to  examine  the  global  identifi- 
ability of  T,j.  Theorem  3.9  uses  the  recursion  relation  of 
Eq  (2.1.14)  to  derive  the  functions  F^,  • • • , 


L 


•jfi. 
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and  then  the  vector  F(t)  is  formed  as  follows: 


'(t)  = [: 


o 1 


n+£.+s- 


J 


(7.2.3) 


By  Theorem  3.9,  if  $ is  open  and  path-connected,  and  if  we 
can  find  a subvector  Fj^(t)  of  F(t]^  formed  from  s components 
of  F(tV  which  is  continuously  differentiable  and  such  that 
3F]^ (t*) /3(()  is  nonsingular  for  some  t*  e(^t^,  tj3  and  all 
e $,  then  ({)  is  globally  identifiable  on  In  this  case 
s = 1 for  4>  = Tj^  and  we  will  consider  ♦ = (0,  »)  so  that  $ 
is  open  and  path-connected.  From  the  recursion  relation 
Eq  (2.1.14) 


^o  = ” ''am(t) 


®’l  = « *am(t) 


^2  = « *am^t> 


(7.2.4) 


X?  ~ ti  am 

n+Us  ■ “ 


3t 


n+£.+s 


Consider  F,  = H x (t).^From  Eg  (7.2.2)  we  obtain 
J-  am 

A (a)T 

''l  - “ I 1)  e + 4j_^(t)]  (7.2.5) 


From  Eqs  (5.4.25)  and  (5.4.26)  we  see  that  e,  (a)  and 

Im 

e (t)  are  independent  of  t^j.  Therefore 


(7.2.6) 
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Using  Eq  (5.4.24)  we  obtain 


n L 


H ='am<‘  I 

n 


Recall  from  Chapter  V that  the  closed  loop  control  matrix 


A (t)  is  defined  by 
c 


A^(t)  = A^(t)  - B3 


A I B(t) 


-X'  ' -1/T, 


where 


X*  = ^2  * * * ^n+nl 


Thus 


3 A^(t) 


This  results  in 


r 0 10' 

fli  = H 1 3X~  ’ 1 


x._(t) 


(7.2.10) 


From  Eq  (7.2.10) , it  is  apparent  that  if 


H * [O 


0 ! l] 


(7.2.11) 


3F, 

then  X— ^ is  nonsingular  (nonzero)  and  continuously  differ- 

entiable  on  the  open  and  path-connected  set  (0,  ®) . There- 
fore we  can  state  that  is  identifiable  on  (0,  <*>)  under 
this  observation.  Also  recall  that 


x^(t)  = 


x(t) 


u(t). 


so  that 


y(t)  = [0  I 1]  x^(t) 


(7.2.11) 


(7.2.12) 


corresponds  to  observing  the  mean  control  U^(t). 

Now  consider  the  simulator  dynamics  described  in  Chapter  VI, 
The  sight  position,  0g,j,»  is  related  to  the  control  by 


the  transfer  function 


K U 
2 RA 


®SA  s(s^  + 3s  + a) 


(7.2.13) 


where  ©sipCt)  = Qs^Ct)  cos0g(t)  with  0 < 0^  < tt/2  and  s > 0. 

From  the  transfer  f\inction,  it  is  apparent  that  there  is  a 
one-to-one  correspondence  between  ©s^iCtQ/  t]  ^RA^^o'  * 

Likewise  there  is  a one-to-one  correspondence  between  the 
functions  E{0g,j,(t)}  and  E{Ujy^(t)}.  Since  is  identifiable 
by  observing  E{Ujy^(t)}  = Ujy^(t)  , then  the  one-to-one  corre- 
spondence implies  is  identifiable  on  (0,  “>)  by  observing 


E{0g^(t)}  = Also  we  know  that  the  target  position 

0^(t)  is  a deterministic  quantity  defined  by  the  predeter- 
mined trajectory.  Therefore  there  is  a one-to-one  corre- 
spondence between  the  functions  0g,j,(t)  and  Gg^(t)  - 0^(t) 
and  also  between  their  mean  values.  This  allows  us  to 
state  that  is  identifiable  on  (0,  <»)  by  observing 
E{0s»p  — 0^}. 

We  will  now  show  that  we  can  come  to  the  same  conclusion 
concerning  identifiability  of  under  the  observations 
E{0st}  and  E{0g,j,  - 0^}  by  continuing  to  examine  the  compo- 
nents F2»  Fjr  and  F4  of  the  recursion  relation.  Using  the 
same  approach  for  F2  = H *x^^  (t)  as  we  used  for  Fj^,  we  obtain 


3F. 


2 _ „ ^ ^am<^> 

3^n  " ^-^n 


A/(t)  = 


3 

Xam(tla)] 

. 3t2  J 

[Ac(t)x3^(t|a)j] 

+ A^^(t)x, 

(7.2.8)  we  can  determine 

A 

- B(t)X'  1 
1 

A B(t)  - i 

- X'A 

^ + ^An  1 

- X B(t)  + 
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‘n 


n 


(7.2.14) 


(7.2.15) 


80  that 


SAj  (t) 


3X' 
B(t)  fE 


n 


3X'  . + 1 m 

t„2 


(l/T^-")B(t) 


If  B(t)  - 2/t^3 
n 


(7.2.16) 


From  Eqs  (7.2.14)  and  (7.2.16)  we  see  that  if  we  observe  a 
mean  state  corresponding  to  a nonzero  component  of 

B(t)r  then  Theorem  3.9  is  satisfied  and  is  observable  on 
(0,  “) . For  the  simulator  described  in  Chapter  VI  we  have 


A^U) 


and 


am 


■ A 

1 B(t)' 

1 

-X' 

1 -lAn 

0 

0 0 0 

1 0 

0 

0 1 C 

1 

1 ° 

0 

-a  -6  C 

' K2  cos0g(t) 

-1 

1 0 C 

1 0 

~^2  "^3 

r • 

E < 

0 

ST 

> 

0 - 0 
ST  T 

^ RA 

(7.2.17) 


(7.2.18) 
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Therefore  on  observation  of  the  form 

y(t)  = [o  0 1 0 o]  x^(t) 


(7.2.19) 


will  result  in  being  observable  on  (0,  ®) . This  corre- 
sponds to  observing  the  mean  sight  acceleration,  E{0  }. 

In  a similar  fashion,  we  find  for  F,  = H x'  „(t)  that 

O cUil 


n n L 


(7.2.20) 


Also 


A^(t) 

c 


A 1 B(t) 


[ A - B(t)X*  ! A B(t)  - l/Tj^B(t) 


L-X'  I L-X'A(t)  + X'  I -X'B(t)  + 1/V 


- A B(t)X'-B(t)X'A+B(t)X'^ 


2 X * A X 

-X'A  + X'B(t)X'  + - :f- 


.2 


4- 


B(t) 


A^  B(t)  - - X'B^(t)  + 


X'B(t)  X'B(t)  1 

-X'A  B(t)  + T_  + T_  *•  _ i 
n n 

n 


A - 
cl  1 

\2 

A 3 1 

A ^ 

L c3  1 

c4  J 

(7,2.21) 


so  that 


-A  B(t)  - B(t)  A I ^ - ||1  B^(t)-^ 

~ “ t_ 
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+ 2 B(t)X'  3X'/8t 
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3A 


c3 


3t 
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3A  - 
c4 

“5? 
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For  the  simulator  of  Chapter  VI 


A B(t)  = 


» 

0 

0 

0 

0 

1 

0 

0 0 10 

0 

0 -a  -3  0 

K2  cos0g(t) 
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rH 
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K2  COS0g(t) 

■3  K2  COS0g(t) 
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(7.2.22) 


(7.2.23) 


From  Eqs  (7.2.20),  (7.2.22),  and  (7.2.23)  we  see  that 
3F2/3Tj^  is  nonzero  and  continuously  differentiable  on  (0,  « 
if  y(t)  =[^01000]  • This  corresponds  to 


observing  the  mean  sight  velocity;  E^g,p).  Therefore  t 


n 


is  identifiable  on  (0,  <*>)  under  this  observation. 


f 

. 

1 


i 

t 

i 


L 


•«* 


For  the  simulator  of  Chapter  VI 


A^BCt)  = 


0 0 0 0 
0 0 10 


t)  -o  -B  0 


-110  0 


COS0g(t) 


-6  Kj  cos0g(t) 


-B  K2  cos0g(t) 

a. 2,21) 

B K2  cos0£(t) 

K2  COS0g(t) 


Eqs  (7.2.24),  (7.2.26),  and  (7.2.27)  indicate  that  an 
observation  of  the  form  y(t)  = [O  0 0 1 O]  will 

result  in  being  identifiable  on  (0,  «»)  . 

In  summary,  we  have  shown  that  Xjj  is  globally  identi- 
fiable on  (0,  «»)  if  the  mean  control  observed 

alone.  Also  for  the  simulator  of  Chapter  VI,  any  of  the 
following  observations  alone  will  also  result  in  x^  being 
identifiable  on  (0,  <») . 

= «i  *am<^^ 


with 


[01000] 

[00100] 

[00010] 


X (t)  = E 
am 


0 - 0„ 

ST  T 


7,2,2  Global  Ident if lability  of  t 


Again  consider  observations  of  the  form 


y(t)  = H x^„(t) 


where 


A^(a)T 


In  this  case  we  will  use  Theorem  3.6  to  examine  the  global 
identifiability  of  t.  That  is,  if  $ is  open  and  path- 
connected  ard  there  exists  a component  of  y(t)  denoted  by 
yj^(t)  which  is  continuously  differentiable  and  3yj^(t*)/3T 
is  nonsingular  (nonzero)  for  t*  e £t^,  1 ' ()>£$, 

then  T is  identifiable  on  $. 

Prom  Eq  (7.2.29) , 


3y(t) 

3t 


3x  (tla)/3T  + A (a)e 


A^(o)t 


(7.2.30) 


Consider  the  second  term  of  Eq  (7.2.30).  From  Chapter  VI 
we  know  the  simulator  dynamics  are  described  by 


A(a)  = 


0 0 0 0 0 

0 0 10  0 
0 -a  -3  0 K2  cos0g(a) 

-1  1 0 0 0 

0000  -1/Xn 


(7.2.31) 


Also 


A.2x2 

e = I + A^T  + — 2^  + 


(7.2.32) 


One  can  inductively  show  that 
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0 
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(7.2.33) 

for  n ^ 2 

Thus  from  Eqs  (7.2.32)  and  (7.2.33)  one  can  state  that  the 


following  elements  of  e 


Ag  (0)t 


will  be  nonzero: 


A (a) T 
.a 


Ell 
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0 

0 
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0 

E21 

^23 
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®25 
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®32 

®33 
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®35 
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E55. 
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Then 


A (a)T 

0 0 0 0 0 

0 0 10  0 

0 -a  -B  0 K2  cosGg 

110  0 0 
0 0 0 0 -1/Tn 

(7.2.35) 

Agj  (a)T 

with  e having  the  form  of  Eg  (7.2.34).  From 

Eg  (7.2.35),  we  see  that  observations  of  the  form 

with 


H = [0  1 0 0 0] 
H = [0  0 1 0 oj 


H = Co  0 0 1 0] 

will  result  in  t being  identifiable  on  (0,  ») . We  cannot 
use  Eg  (7.2.35)  to  state  identif lability  under  the 
observation  H=Co  0 0 0 l]  because  the  mean  estimator 

error  of  the  control  (i.e.,  the  fifth  component  of  ©^^^^^(a)) 

\ 
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is  zero.  However  we  can  show  identif lability  under  this 
observation  by  referring  to  the  transfer  function  for  the 
simulator. 

Recall  that  the  sight  position,  is  related 

to  the  control  U by 
RA 

K U 

®RA  = (7.2.36) 

- s(s^  + Bs  + a) 

where  0g^(t)  = 9g^(t)  cos0g(t)  with  0 < 0^  < ir/2  and  s > 0. 

Thus  there  is  a one-to-one  correspondence  between  the 

functions  0g^(t)  and  and  also  between  the  ensemble 

averages  of  these  functions;  i.e,,  E{0  (t)}  and  E{U  ,(t)}. 

ST  RA 

since  0^(t)  is  a deterministic  function  independent  of  t, 

this  allows  us  to  state  that  there  is  a one-to-one  corre- 
.spondence  between  the  functions  E{0  - 0 } and  E{u  (t)}. 


Therefore,  since  x is  identifiable  on  (0,  «»)  by  observing 


E{0gq,(t)  - 0.p(t)},  we  can  state  that  x is  identifiable  on 
(0,  «»)  by  obseirving  E{U  (t)}. 

In  summary,  x is  identifiable  on  (0,  ®)  for  the  simu- 
lator system  described  in  Chapter  VI  under  any  of  the 
following  observations 

where 

% = Co  1 0 0 o] 

= [o  0 1 0 O] 

H3  = [O  0 0 1 0] 

H4  = [0  0 0 0 1] 


7.2.3  Global  Identification  of 

Let  P (t)  be  the  variances  of  the  states  x, (t) . Now 
axx  a 

consider  observations  of  the  form 

y(t)  (7.2.37) 

To  excimine  the  global  f?Jentif lability  of  p^,  the  proportion- 
ality constants  for  calculating  observation  noise  variance 
kernels,  we  will  again  use  Theorem  3.9.  Let  $ be  the  open 
and  path-connected  set  (0,  <»)  and  again  let 

F(t)  - [f/  Fi’'  • ■ • (7.2.38) 
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where 


Fq  = y(t) 
Fj  = y(t) 
F,  = y(t) 


.n+A+s 
3 y(t) 


n+£+s  " 

p£  is  identifiable  on  4)  if  we  can  find  a component  of  F(t), 
say  Fj^Ct),  which  is  continuously  differentiable  and 
3Fj^(t*)/3p^  is  nonzero  for  some  t*  e t^j  and  all 

Pi  ^ (0,“) . From  Eq  (5.4.32),  we  know 


A^(o)t  A^'^(a)T  A ’’(a)T 

+ P^e  ® 


P^(t)  = e ^ ?2 (a)e  ^ 


Aa(cr)T  rj. 


+ e P3  (t)  + P4(t)  + P5(t) 


(7.2.39) 


Using  Eqs  (5.4.28)  to  (5.4.31)  we  obtain  the  following  as  a 
first  step  to  evaluating  8Fj^  (t) /9p . 


3Pa(t)  g I"  Ag(o)T.  Ag'^(a)T  . A^'®’(a)T 
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(t)  _ 9 r 
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(7.2.40) 


(7.2.41) 
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(7.2.42) 
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(7.2.44) 


T 3K  T 3V  ^ 

^ s + e Pi  C_  aT C P e 

i 3p.  1 a 3P^  a 1 


(7.2.45) 


Now  consider  the  simulator  of  Chapter  VI  where 


C = 
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0 0 0 1 0 
-110  0 0 


Thus 


T -1 
C V C 
a a 
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■110  0 0 


Now  which  is  a term  of  ^2'  will  have  the 

form  PiC/v‘^C^Pl  - 


®11^2 

®12^2 

®13'^2 

1-1  2-1 
®14^1  ■^^14'^2 
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From  Eq  (7.2.34)  we  know 
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43 
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^55 

A ^a"*^  T -1  T ^a*^^ 

Now  consider  J = e PiC,  V C,P,  e . From  the  above, 

1 ^a  a 1 ' 

it  is  seen  that  the  diagonal  element  will  have  the  form 


-r  -1  . 2„  -1 

J = e V + e V 
44  44  1 44  2 


* J-  V .J 


with  e 


Combining  the  above  results  with  Eqs  (7.2,40)  through 
(7.2.45)  we  see  that  aPa(t)/3pi  will  have  a nonzero 
diagonal  term  3Pjj44  * "^his  means  observations  of 

the  form 


with  H 


will  result 


will  be  nonzero  for 


Therefore,  for  the  simulator  dynamics  described  in 
Chapter  VI,  the  parameters  p.  (i  = 1,  2)  are  identifial 


on  the  open  interval  (0,  <»>)  under  an  observation 


This  implies  that  3^a44  nonzero.  Also  3P55(t)/3Pj, 

is  nonzero.  Using  Theorem  3.9  and  the  same  approach  as  for 
we  conclude  that  is  identifiable  on  (0,  »)  under  the 
observation 


y(t) 


with  H = [0  0 0 1 O]. 


Again  this  corresponds  to  observing  the  standard  deviation 
of  the  tracking  error. 

Also  Pjn  is  identifiable  on  (0,  <»)  under  the  observation 
y(t)  =. 

with  H = [O  0 0 0 l]. 

This  corresponds  to  observing  the  standard  deviation  of 

•• 

the  sight  acceleration,  ©cm,  or  the  standard  deviation  of 


the  control,  Up^. 

Again  recall  that  the  sight  position  0g^  is  related  to 
the  control  Up,  by  the  transfer  function 


^2  V 

, 2 

s (s  + 3s  + a) 


where  ©s^^^)  “ ©g^  (t)  cosG^  (t ) with  0 < 0^  < tt/2.  This 

implies  a one-to-one  correspondence  between  the  functions 

0 (t)  and  U (t) . Also  this  assures  a one-to-one  corre- 

ST  RA 

2 2 

spondence  between  the  functions  0g,p(t)  and  Ujy^(t)  and  thus 


between  the  standard  deviation  of  0^  and  U . 

ST  RA 


Now  since 


Pjj^  is  identifiable  under  the  observation  of  the  standard 


deviation  of  (t) , we  also  can  state  that  p is  identi- 
RA  in 

fiable  under  the  observation  of  the  standard  deviation  of 

Since  the  target  position  is  a deterministic  quan~ 

tity,  the  standard  deviation  of  the  tracking  error, 

©_(t)-0  (t),  is  equal  to  the  standard  deviation  of  the 
ST  T 

sight  position.  Thus  is  also  identifiable  under  the 
observation  of  the  standard  deviation  of  the  tracking 
error. 

Using  a similar  argument,  we  can  state  that  is 

identifiable  on  (0,  «)  by  observing  the  standard  deviation 

of  the  function  (t) . 

RA 

In  summary,  both  p^(i  = 1,  2)  and  Wj^  are  identifiable 
on  (0,  “)  by  observing  either  the  standard  deviation  of 

V“>  ®ST*‘’  - 

7.2.5  Global  Identification  of 

Recall  that  in  Section  5.4  an  alternative  model  for 
the  plant  disturbance  was  presented.  The  alternative  has 
the  white  noise  component  equal  to  zero  and  thus  the  plant 
disturbance  is 


w(t)  = W2<t)  (7.2.47) 

and  the  variance  of  the  disturbance  is  modeled  by 

W(t)  = 2 W2^(t)  (7.2.48) 

where  may  be  viewed  as  a correlation  time  for  W2(t). 

One  could  assign  a value  to  equal  to  an  estimate  of  the 
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correlation  time  or  one  could  attempt  to  estimate  its  value 
using  experimental  data.  In  case  of  the  latter,  we  are 
Interested  in  this  pareimeter's  identifiabilit^'. 

Following  the  Scime  approach  as  in  preceding  sections 
we  find 


where 


^^1 ^ „ 3W(t)  T 

3Tc  a 3t^  a 


3W(t) 
3(t)  " 


2 

2 w,  (t)  I 

_ I 

I 

0 I 


This  term  couples  into  term 


Af(t)  = A^(t)  ~ PiCg'^V  ^(a) 


(7.2.49) 


f 


Using  the  same  approach  as  for  and  Wj^  we  can  show  that 

1 is  identifiable  on  (0,  “)  under  an  observation  of  the 
c 

standard  deviation  of  the  tracking  error  0g^  - 0^  or  of 
the  control,  Uj^. 


7.2.6  Svunmary  of  Global  Identifiability 

In  the  previous  sections  we  examined  the  identifiability 
of  the  optimal  control  model  parameters  under  observations 
of  the  form 
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and  standard  deviation  of  these  states. 

Foi  the  simulator  described  in  Chapter  VI,  we  show  that 
the  parameters  were  identifiable  on  (0,  «)  under  the  obser- 
vations indicated  by  the  following  matrix: 


Mean 

States 

standard 

Deviations 

• 

0 

ST 

• • 

0 

ST 

0 -0 
ST  T 

V 

0 -0 
ST  T 
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Tn 

X 

X 

X 

X 

T 

X 

X 

X 

X 

Pi 

X 

X 

P2 

X 

X 

Wi(Tc) 

X 

X 

Pm 

X 

X 

Thus  from  the  theory  developed  in  Chapter  III,  we  know  that 
4>  = t P2  identifiable  under  obser- 

vations of  the  form  of  Eq  (7.2,50)  with 

H = C0  0010|  00010]  or  H = [0  000lj  0000  l]. 

7.3  Computational  Results 

In  the  preceding  sections  we  have  shown  that  the  param- 
eters of  the  optimal  control  model  can  be  identified  from 
observations  of  the  mean  and  standard  deviation  of  either 
the  tracking  errors  or  the  manual  control  inputs  to  the 
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tracking  system.  In  this  section  the  results  of  compu- 
tations which  were  performed  to  identify  the  parameter 
values  will  be  presented.  We  obtain  data  from  actual 
system  operation  using  the  simulator  system  described  in 
Chapter  VI.  This  yields  the  data  shown  in  Figs.  6.8.1 
through  6.8.4,  which  results  from  ensemble  averaging  many 
time  histories  of  tracking  errors.  This  simulator-gener- 
ated data,  denoted  by  z(t^),  t^j^  e [t^,  t^]  , has  corre- 
sponding data  resulting  from  model  calculations  which  are 
denoted  by  y (t^^)  . We  want  to  find  an  estimate  of  the  param- 
eter  vector  c|)  which  minimizes  the  cost  functional 

k qv 

[z(t.)  - y(t.)]  [z(t.)  - y(t.)]  (7.3.1) 

i=l  ^ ^ 

T 

where  <|)  = [t„  t P2  ^m  * 

The  measured  data  z(t^)  to  be  used  in  the  identification 
computations  consists  of  the  mean  and  standard  deviations 
of  the  tracking  errors.  As  discussed  in  the  previous 
sections,  this  is  sufficient  to  allow  identification  of  the 
parameter  vector  <|> . Other  sets  of  measured  data  could  be 
considered,  such  as  the  mean  and  standard  deviation  of  the 
control  input  or  combinations  of  tracking  error  and  control 
input  data.  In  this  research,  we  will  confine  our  compu- 
tations to  those  using  tracking  error  data;  follow  on 
research  could  make  use  of  other  possibilities  for  measured 
data  sets. 
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We  will  use  the  gradient  method  for  minimizing  the 

cost  functional  as  discussed  in  Section  4.4.1.  This 

requires  calculation  of  the  gradient  equations  9y(t)/3<|). 

Since  the  measured  data  is  a function  of  the  mean  states, 

X _(t),  and  their  covariances,  P, (t) , we  need  to  compute 
am  a 


am 

d(|> 

and 


9x  (t)  I 
am  I 

9(|),  j 


9x  (t)  I 
am  I 

Ho  ! 


I 9x 
I am 


9(]> 


6 J 


(7.3.2) 


3P^ (t) 

cl 

9^ 


3P  (t) 
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9(|)- 


I 9P  (t)  I 

I a __  I 

! H.,  ! 


I 3P  (t) 
I a 

: H. 


(7.3.3) 


This  is  a very  tedious  operation  that  follows  in  a straight- 
forward manner  from  the  equations  for  propagating  the  mean 
states  and  covariances.  Appendix  F summarizes  the  results 
of  gradient  equation  evaluation  for  the  simulator  under 
discussion. 

Computer  programs  were  developed  to  propagate  the  mean 
states  and  covariances  described  by  the  coupled  matrix 
differential  equations  in  Section  6.5.  Also  the  gradient 
equations  developed  in  Appendix  F are  solved  using  the 
digital  computer.  Appendix  G presents  a flow  diagram  of 
the  programs  and  subroutines  used  to  make  the  required  com- 
putations . 

The  gradient  method  was  used  to  adjust  the  parameter 
vector  (|)  treating  each  axis  as  an  independent  case.  Recall 
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that  the  particular  simulator  used  in  this  research  has  two 
people  in  the  tracking  system;  one  operates  the  azimuth 
control  and  one  operates  the  elevation  control.  This  two- 
operator  system  supports  the  assumption  of  treating  the  axes 
separately  (other  than  the  cos0^(t)  term  appearing  in  the 
traverse  equations) ; however,  it  is  likely  that  some 
coupling  effects  do  occur.  This  will  be  discussed  in  con- 
junction with  some  of  the  results  to  be  presented  and  could 
be  addressed  in  future  research  as  an  important  refinement 
to  the  model. 

In  the  next  sections  the  results  obtained  from  the 
parameter  identification  computations  will  be  presented  as 
well  as  a comparison  of  model  and  simulator  data.  Recall 
that  two  trajectories  were  considered;  the  profiles  are 
shown  in  Figs.  6.7.1  and  6.7.2  respectively.  As  discussed 
in  Section  6.7,  the  acquisition  phase  takes  place  during 
the  first  10  seconds  of  the  trajectory.  This  portion  of 
the  trajectory  was  avoided  for  identification  purposes 
because  the  assumption  of  the  least  squared  error  cost 
functional  for  the  trackers  would  be  suspect.  To  allow  any 
transients  remaining  from  the  acquisition  phase  to  damp  out, 
the  identification  calculations  were  started  at  15  seconds 
into  the  trajectories  and  were  computed  over  10  seconds  to 
a final  time  of  25  seconds.  This  covered  the  crossover 
condition  (x  = 0)  for  both  trajectories  and  was  considered 
adequate  to  provide  a basis  for  comparing  model  and 
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simulator  data,  as  well  as  for  the  gradient  pareuneter 
identification  procedure. 


I 


7.3.1  Trajectory  1,  Azimuth  Axis 

An  initial  value  for  the  vector  was  assumed  which 
corresponds  to  nominal  values  used  by  Kleinman,  et  al 
(Refs  5r  20,  21,  22,  23).  These  values  were 

Tn  = .1 

T = .2 

Pj^  = .0314 

p2  = .0314 

P™  = -01 
m 

= .000001 

Note  that  the  plant  disturbance  noise,  Wj^,  was  introduced 
as  a new  parameter  in  this  research  to  account  for  system 
noise  and  model  uncertainities , so  that  no  nominal  value 
was  available  from  the  literature.  The  initial  estimate 
for  corresponded  to  an  estimate  for  a nominal  value  of 
1^0^  (t)J^  for  the  trajectory. 

Recall  from  Eq  (4.4.10)  that  the  change  in  (P  for  each 
iteration  is  given  by 

lc+1  k f kl-1 

A<|>  = <l>  - [n  J 


3J(j)^) 


(7.3.4) 


where  N is  a symmetric  positive  definite  matrix  to  adjust 
the  step  size  with  the  following  form 


1R2 


0 


0 


2 

0 Nj’' 

0 0 0 
4 

0 0 0 N. 

5 

0 0 0 0 


= I [z(t.)  - y(t.)]  [z(t.)  - y(t.)]  (7.3.5) 

i=l  ^ 1 1 1 


I Lz(t.)  - y(t.)]  (7.3.6) 

i=l  9<|)  J 


As  stated  above  z(t^)  is  a two  element  vector 

'zi(ti) 

z(t  ) = 

Z2(ti) 

where 

Zi(ti)  = measured  mean  tracking  error  at  t^ 

= standard  deviation  of  measured  tracking  error 
at  tj^ 

Table  7.3.1  summarizes  the  results  of  using  the  gradient 
technique  to  adjust  the  values  of  <t>  for  the  azimuth  axis  of 
trajectory  1. 

The  following  observations  are  made  from  Table  7.3.1: 
(1)  The  cost  functional  J((l)  ) is  reduced  at  each 


(7.3.7) 


Table  7.3.1 


r 


SUMMARY  OF  GRADIENT  ITERATIONS 
TRAJECTORY  1,  AZIMUTH  AXIS 


k=l 

1 3J(<>^) 

2 3,j,k 

2[n*']  ^ 

J((J)^)  = . 946x10"^ 

.100 

.1031x10"^ 

50 

|m((|)^)  | = . 579x10“^® 

.2000 

.5013x10"? 

50 

.0314 

.4082x10";: 

50 

.0314 

.9632x10“^ 

50 

.0100 

.2492x107^ 

50 

.000001  .2009x10^  5x10“® 


k=2 


J(4>^)  = . 335x10“'^ 

.105 

.4236x10"? 

50 

|m(<|»^)  | = .469x10"2° 

.2025 

.3307x10“^ 

50 

.0515 

.6250x10", 

50 

.0362 

.1122x10":: 

50 

.0225 

.1017x10", 

.00001104 

.2084x10"'^ 

5x10"® 

k=3 

J(<|»^)  = .303x10“4 

.107 

.1028x10"? 

10 

|m(4)^)  | = . 375x10"^° 

.2042 

.0546 

.1367x10"2 

.3208x10"^ 

10 

10 

.0418 

.4888x10“? 

10 

.0283 

.3034x10"^ 

.00001105 

-.1019 

10  ® 

k=4 

J (4)^)  = .288x10"^ 

.108 

.1736xl0“f 

|m((|)4)  | = . 339x10"^° 

.2055 

.7725x10“® 

.0578 

.2469x10“? 

.0467 

.3088x10"; 

.0313 

.1478x10"^ 

.00001003 

-.1404 

I 


iteration  indicating  the  gradient  method  for  adjusting  41 
was  functioning  properly. 


-1 


C2)  The  values  of  the  diagonal  of  [n'^]  were  selected 
off  line  to  give  reasonable  step  sizes.  Since  the  cost 


functional  was  reduced  at  each  step,  this  indicates  the 


choices  for  [n^]  were  satisfactory.  Recall  from 


Section  4.4  that  too  large  a step  size  can  cause  overshoot 
and  an  increase  in  the  cost  functional. 

(3)  The  values  of  the  gradients  9J/9(|)  decreased  with 
each  iteration  as  one  expects.  It  v/as  found  that  after 
four  iterations  the  gradients  were  sufficiently  small  to 
make  the  change  in  <j)  very  slight.  The  stopping  condition 
was  based  on  | | A(|)^  | | being  less  than  a preselected  small 
value . 


(4)  Recall  from  Eq  (7.1.1  that  the  parameter  vector  4> 

N T 

is  locally  identifiable  if  M(4>)  = 1 C 9y  (t . ) /9(})] 

2 _ 1 1 


[3y(tj^)/H]  is  nonsingular. 


i=l 


k k 

The  determinant  of  M(<{)  ),  denoted  by  |m(<J)  )|  was  evaluated 


at  each  iteration  and  appears  in  Table  7.3.1.  In  each 

Lk 


case  |m(<J)  ) I is  nonzero*  indicating  that  M((J))  is  nonsingu- 
lar. This  supports  the  hypothesis  that  the  parameter 


vector  (fi  is  locally  identifiable  at  each  of  the  iteration 

ik 


values,  <t>  . This  suggests  the  possibility  of  using  a 
second  mode  of  Newton  Raphson  iterations  to  converge 
rapidly  near  the  end. 

Fig.  7.3.1  compares  the  model  data  to  that  obtained  from 
the  simulator  for  the  parameters  from  the  fourth  Iteration 
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of  the  azimuth  axis  identification  for  trajectory  1.  The 
model  mean  tracking  error  is  slightly  less  than  the  simu- 
lator error  prior  to  target  crossover  at  20  seconds.  How- 
ever, the  model  results  agree  with  the  simulator  in  that 
the  sight  tends  to  lead  the  target  prior  to  crossover.  This 
supports  the  configuration  of  the  model  treating  target 
angular  acceleration  as  a plant  disturbance  as  opposed  to 
the  third  derivative  of  angular  position  for  example.  Under 
this  configuration,  the  operator  (as  modeled  by  the  Kalman 
filter/predictor)  estimates  the  target  velocity  based  on 
the  assumption  that  acceleration  is  a zero  mean  random 
disturbance.  Fig.  7.3.2  shows  the  actual  target  angular 
velocity  and  acceleration  as  a function  of  time. 

The  value  of  (strength  of  the  white  noise  portion 
of  the  plant  disturbance)  accounts  for  system  noise,  model 
uncertainities , and  the  ability  of  the  operator  to  adapt  to 
different  trajectories.  Increasing  will  result  in  the 
Kalman  filter  gain  increasing  and  thereby  tend  to  increase 
the  weight  of  current  observations  of  the  operator.  Thus 
the  value  of  together  with  the  operator's  estimate  of 
the  value  of  the  angular  acceleration,  reflect  the  adapt- 
ability of  the  model  filter  and,  as  such,  the  operator 
adaptability.  Such  techniques  have  been  used  in  other 
applications  as  discussed  in  Refs  15,  16,  and  34. 


■( 


, ■ 


■mr 
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Table  7.3.2  presents  a stairanary  of  the  gradient  calcu- 
lations for  the  elevation  axis  of  trajectory  1.  Again  the 
cost  functional  is  reduced  at  each  iteration.  Moreover 
|m(<})  ) I is  nonzero  indicating  local  identif lability . 

Figure  7.3.3  shows  the  model  and  simulator  response  with 
the  model  data  corresponding  to  the  parameter  values  after 
the  third  iteration  of  the  identification  procedure.  From 
the  figure  it  is  apparent  that  the  elevation  tracking  task 
is  not  very  difficult  for  trajectory  1.  The  mean  response 
of  the  simulator  and  the  model  agree  in  that  the  sight  tends 
to  be  above  the  target.  Again  this  tends  to  support  the 
model  configuration  with  target  acceleration  as  a plant 
disturbance.  The  actual  target  angular  velocity  and  accel- 
eration for  the  elevation  axis  of  trajectory  1 is  shown  in 
Fig.  7.3.4.  The  mean  tracking  error  introduced  by  the  tar- 
get acceleration  in  the  model  response  results  in  general 
agreement  of  model  and  simulator  response. 

The  simulator  mean  tracking  error  has  a relatively 
small  oscillatory  characteristic  near  crossover  not  present 
in  the  model  response.  This  could  be  attributed  to  the 
changing  image  of  the  target  on  the  tracking  monitor  during 
the  crossover  phase,  with  the  operator  adjusting  the  ele- 
vation cross  hair  to  account  for  what  is  perceived  as  a 
small  change  in  elevation.  For  a deterministic  trajectory, 
this  would  be  a consistent  effect  for  a given  trajectory 
and  therefore  appear  as  a perturbation  to  the  mean  tracking 


Table  7.3.2 


SUMMARY  OF  GRADIENT  ITERATIONS 
TRAJECTORY  1,  ELEVATION  AXIS 


k=l 

1 3J(<|)^) 
^ 3(|)^ 

2[n>'] 

J((J)^)  = . 461x10“^ 

.100 

.8140x10"^ 

50 

|M(4,^)  | = . 465x10"^^ 

.2000 

.1064x10"^ 

50 

.0314 

.1034x10"3 

50 

.0314 

.4654x10"'^ 

50 

.0100 

.2331x10“3 

50  c 

.00001 

.6021 

.5x10"^ 

k=2 

J(({)^)  = . 378x10“^ 

.104 

-.5019x10"^ 

10 

|M((j,^)  1 = .292x10"2® 

.2005 

-.5228x10"^ 

10 

.0366 

-.5203xl0"f 

10 

.0337 

-.1139x10"^ 

10 

.0216 

-.3911x10“^ 

10 

.00001301 

-.5538 

.5x10“5 

k=3 

J(4>^)  = . 370x10“^ 
|m((|)^)  | = . 290x10"^® 

.104 

.2000 

-.5425x10“^ 

-.5119x10“^ 

.0366 

.0327 

.0214 

.00001352 

-.4673x10"^ 

-.1144xl0“f 

-.4150x10"^ 

.5293 

EL  ERROR  (mil) 

30.00  -15.00  0.00  15.00  30.00 


error  as  indicated  by  the  simulator  response  in  Pig.  7.3.3. 
During  this  phase  of  the  trajectory  the  average  position 
of  the  azimuth  sight  is  also  changing  from  a sight  leading 
situation  to  a sight  lagging  situation.  The  position  of 
the  azimuth  sight  on  the  finite  image  size  may  have  some 
effect  on  the  elevation  operator  positioning  his  cross  hair. 

7.3.3  Trajectory  2,  Azimuth  Axis 

As  seen  in  Fig.  6.7.2,  trajectory  2 requires  the  tar- 
get to  maneuver  more  than  the  straight  and  level  fly-by  of 
trajectory  1.  This  is  also  seen  by  the  velocity  and  accel- 
eration time  functions,  which  are  shown  in  Fig.  7.3.5  for 
trajectory  2,  azimuth  axis.  Table  7.3.3  presents  a summary 
of  the  gradient  calculation  results.  Again  the  cost  func- 
tional is  reduced  at  each  iteration  and  the  determinant  of 

k k 

M((|)  ) is  nonzero  at  each  iteration  value  of  ({>  . Iteration 

one  used  the  values  of  (p  arrived  at  for  the  last  iteration 

of  the  azimuth  axis  for  trajectory  1.  The  most  significant 

change  in  the  adjustment  of  (P  during  the  four  iterations  on 

the  trajectory  was  to  increase  the  strength  of  the  white 

noise  portion  of  the  disturbance,  This  supports  our 

previous  comment  that  reflects  the  adaptive  capability 

of  the  operator.  An  increase  in  will  increase  the  gain 

of  the  Kalman  filter  in  the  model  and  result  in  the  filter 

increasing  the  weight  of  recent  data.  This  is  desirable 

when  target  is  maneuvering  significantly,  as  is  the  case 

for  this  trajectory.  Of  course,  some  of  this  adaptability 
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1 f Table  7.3.3 


SUMMARY  OF  GRADIENT  ITERATIONS 
TRAJECTORY  2,  AZIMUTH  AXIS 


k«l 

1 3J(d)^) 

2 

J((>^)  = .277x10"3 

|=. 189x10”^^ 

.108 

.2055 

.0570 

.0467 

.0313 

.00001003 

.2156x10"3 

-.2217x10”3 

.4254x10“3 

-.1566x10“3 

.8054x107^ 

.1151x10-^ 

k=2 

J(4>^)  = . 236x10"^ 

|M((j)2)  | = .896x10"^3 

k=3 

.110 

.2033 

.0621 

.0451 

.0394 

.00002154 

“. 1449x10"! 
-.5845xl0";r 
.1678x10", 
-.9019x10"^ 
-.7995x10"^ 
.9472 

J(<J>^)  = .201x10"3 

|m((J)^)  I = . 486x10"^^ 

k=4 

.109 

.1975 

.0623 

.0361 

.0386 

.00003101 

-.2627x10"^ 

-.3918x10"; 

.9452x10"^ 

-.6307x10"; 

.8722x10"^ 

.7853 

J(4»^)  = . 199x10“^ 
|M(«j)^)  | = . 471x10"^^ 

.109 

.1971 

.0624 

.0355 

.0387 

.00003180 

-.2367x10"!? 

-.3848x10"; 

.9431x10", 

-.6225x10"^ 

.8641x10“^ 

.7740 

reflected  by  the  model  also  having  the  square  of  the  angu- 

2 

lar  acceleration,  » contribute  to  the  strength  of 

the  disturbance  as  seen  by  the  filter. 

Fig.  7.3.6  compares  the  model  response,  using  parameter 
values  after  iteration  four,  to  the  simulator  response. 

Here  the  mean  tracking  error  of  the  model  is  of  slightly 
greater  magnitude  than  that  of  the  simulator,  although  both 
have  the  sight  leading  the  target.  The  larger  standard 
deviation  of  the  simulator  between  15  and  18  seconds  may  be 
attributed  to  crossover  effects  from  the  elevation  axis. 

As  seen  from  Fig.  7.3.7,  the  elevaci  acceleration  during 
this  period  is  relatively  high,  which  . ds  to  increase  the 
variance  of  the  elevation  tracking  error.  The  azimuth 
operator  may  be  affected  by  this  apparent  disturbance  in 
the  elevation  axis  with  the  result  being  an  increase  in  the 
azimuth  tracking  error  variance.  One  could  have  the  model 
account  for  this  by  one  of  two  methods:  (1)  add  a term  to 

equivalent  observation  noise  in  the  azimuth  axis  that  is 
proportional  to  variance  of  the  observed  quantities  in  the 
elevation  axis  (and  vice-versa)  or  (2)  add  a term  in  the 
variance  of  the  plant  disturbance  seen  by  the  Kalman  filter 
in  the  azimuth  axis  that  is  proportional  to  the  square  of 
the  elevation  axis  angular  acceleration  (and  vice  versa) . 
Either  of  these  modifications  would  provide  an  increase  in 
the  variance  of  the  tracking  error  that  was  proportional  to 
the  variance  in  the  other  axis. 
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FICRJRE  7.3.6 

TRAJECTORY  2,  AZIMUTH  AXIS 
MODEL  AND  SIMULATOR  TRACKING  ERRORS 


7.3.4  Trajectory  2,  Elevation  Axis 


The  elevation  tracking  task  is  also  more  difficult  for 
trajectory  2 than  for  trajectory  1.  This  is  exemplified  by 
the  angular  velocity  and  acceleration  time  functions  shown 
in  Fig.  7.3.7.  Table  7.3.4  summarizes  the  gradient  calcu- 
lations for  this  trajectory  and,  as  with  the  other  cases, 
the  cost  functional  is  reduced  at  each  iteration  and  the 
determinant  of  M(<|)  ) is  nonzero,  supporting  local  identifi- 
ability  of  4»  . Fig.  7.3.8  compares  the  model  response  to 
the  simulator  response  using  the  parameter  values  after 
iteration  four.  The  model  mean  and  standard  deviation  are 
generally  in  good  agreement  with  the  corresponding  simu- 
lator values.  The  simulator  mean  tracking  error  is  less 
accurate  than  that  predicted  by  the  model  in  the  19  to  23 
second  region  when  the  sight  is  changing  from  an  average 
leading  position  to  one  in  which  the  sight  lags  the  target. 
From  Fig.  7.3.7,  it  is  seen  that  this  is  the  region  in  which 
the  angular  acceleration  is  changing  from,  a negative  to  a 
positive  value.  The  overshoot  in  the  22  second  region 
exhibited  by  the  simulator  could  be  attributed  to  the  oper- 
ator anticipating  the  point  in  the  trajectory  where  the  tar- 
get angular  velocity  is  no  longer  increasing  in  magnitude. 
Although  the  operators  were  presented  with  one  of  three 
trajectories  at  random,  it  is  suspected  that  after  several 
trials  the  characteristics  of  the  trajectories  are  learned 
to  some  degree.  Another  possible  explanation  for  the  dif- 
ference is  that  the  numerically  calculated  acceleration 
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Table  7.3.4 


SUMMARY  OF  GRADIENT  ITERATIONS 
TRAJECTORY  2,  ELEVATION  AXIS 


k=l 

1 3J((|)^) 
^ 3(fr^ 

2 [n*^J 

J((|)^)  = . 211x10“^ 

.104 

.4029x10"3 

10 

|m((|)^)  | = .381x10"15 

.2000 

.2773xl0~i 

10 

.0366 

.3984x10"; 

10 

.0327 

.7430x10"; 

10 

.0214 

.9502x10"^ 

10  . 

.00001352 

.4070 

10"5 

k=2 

J(4>^)  = . 136x10"^ 

.108 

.2816x10"^ 

10 

|M(4>^)  | = .435x10“^5 

.2028 

.0406 

.2187x10"3 

.2792xl0”f 

10 

10 

.0401 

.4623x10"^ 

10 

.0309 

.6433x10"3 

10  . 

.00001752 

.2070 

10"5 

k=3 

J((fr3)  = . 104x10"^ 

.111 

.1816x10"3 

10 

|m((|)^)  1 = . 487x10"^^ 

.2049 

.0434 

.1535x10"3 

.1842x10"3 

10 

10 

.0447 

.2929x10"3 

10 

.0373 

.3928x10"3 

10 

.00001959 

.1243 

10"5 

k=4 

J((|)^)  = .910x10"4 

.113 

.1148x10"3 

|m(4»^)  1 = . 518x10"^^ 

.2064 

41452 

.1021x10"; 

.1178x10", 

.0476 

.1877x10"^ 

.0412 

.2227x10"3 

.00002083 

.8173x10-1 
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FIGURE  7.3.8 

TRAJECTORY  2,  ELEVATION  AXIS 
MODEL  AND  SIMULATOR  TRACKING  ERRORS 


la 

- Model  Mean 


X X X Simulator  Mean 


disturbance  may  not  exactly  duplicate  the  actual  angular 
acceleration  seen  by  the  simulator.  The  possibility  of 
some  crosscoupling  between  the  two  axes  could  also  be  a 
factor.  It  is  noted  that  during  the  19  to  23  seconds,  both 
the  azimuth  and  elevation  operators  have  relatively  diffi- 
cult tracking  tasks.  The  effect  of  a changing  sight 
position  in  azimuth  axis  on  the  elevation  operator,  and 
vice  versa,  is  not  readily  apparent;  however,  some  cross- 
coupling may  take  place.  As  mentioned  before,  some  small 
effects  could  be  introduced  by  the  changing  target  image  on 
the  TV  monitor. 

7,3.5  Alternative  Model  for  Plant  Disturbance 

As  discussed  in  Section  5.4,  one  could  model  the  plant 
disturbance  without  the  white  noise  component;  that  is,  the 
disturbance  would  be  only  the  deterministic  quantity  W2 (t) . 
In  our  case  W2 (t)  is  the  angular  acceleration.  Following 
the  approach  by  Kleinman  (Refs  20,23),  the  covariance  kernel 
of  this  disturbance  as  seen  by  the  Kalman  filter  is  modeled 
by 

E{[w(tj^)  - W(tj)]'^w{t2)  - w(t2)]}  = 2 0(tj)6(t^  - t2) 

where  x may  be  viewed  as  a correlation  time  for  0(t). 
c 

A second  set  of  calculations  was  performed  v’ith  this 
alternative  model  for  the  disturbance,  x^  replaced  Wj^  as 
a quantity  to  be  identified  in  the  parameter  vector  <j). 

Table  7.3.5  and  Figs.  7.3.9  through  7,3,12  summarize  the 
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Table  7.3.5 


RESULTS  OF  GRADIENT  CALCULATIONS 
ALTERNATIVE  MODEL  FOR  PLANT  DISTURBANCE 


Trajectory 

Final 

Iteration 

J (<!)’') 

_ 1 

1, 

Azimuth 

4 

.119 

.291x10"^ 

-.2798x10"^ 

.2078 

-.2590x10"5 

.0977 

-.4452xl0"f 

.0572 

-.1637x10"5 

.0579 

-.1584xlO~Z 

6.919 

-.1417x10"' 

1, 

Elevation 

3 

.110 

.309x10"^ 

-.5943x10"^ 

.2044 

-.5370x10"’ 

.0587 

-.3171x10"^ 

.0503 

-.8055xl0"f 

.0416 

-.8839x10"^ 

4.234 

-.1212x10"^ 

2f 

Azimuth 

4 

.109 

.162x10"^ 

-.1410x10"^ 

.2026 

-.2552x10"^ 

.0551 

.5174x10"^ 

.0466 

-.3243x10"^ 

.0385 

1.303 

.1588xl0"J 

.6521x10"'’ 

2, 

Elevation 

4 

.114 

.120x10"3 

.1004x10"^ 

.2061 

.1017x10"3 

.0602 

.7093x10"^ 

.0531 

.1599x10"^ 

.0453 

.7540x10"4 

5.099 

-.1672x10"5 

I 
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FIGURE  7.3.11 
TRAJECTORY  2,  AZIMUTH 
MODEL  AND  SIMULATOR  TRACKING  ERRORS 
(ALTERNATIVE  MODEL  FOR  PLANT  NOISE) 


X 


Model  Mean 


Simulator  Mean 


FIGURE  7.3.12 
TRAJECTORY  2,  ELEVATION 
MODEL  AND  SIMULATOR  TRACKING  ERRORS 
(ALTERNATIVE  MODEL  FOR  PLANT  NOISE) 
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results  obtained  using  this  approach.  Comparing  this  data 

to  that  obtained  with  the  first  model  for  plant  noise,  we 

see  that  the  cost  functional  obtained  after  the  same  number 

of  iterations  is  not  significantly  different  in  either  case. 

A somewhat  better  agreement  of  model  and  simulator  data  is 

obtained  for  the  new  formulation  of  the  plant  disturbance 

for  the  elevation  axis  of  trajectory  1 and  the  azimuth  axis 

of  trajectory  2,  although  both  provide  acceptable  results. 

A model  could  be  formulated  which  includes  the  white 

noise  components  Wj^(t)  in  the  disturbance  and  also  maintains 

the  factor  t in  the  variance  kernel.  That  is 
c 

..2 

W(t)  = Wj^(t)  + 0 (t) 

and  the  process  noise  covariance  kernel  as  seen  by  the 
Kalman  filter  is 

..2 

W(t)6(t)  = Wj^fi(t)  + 2 0 (t)6(t) 

In  this  manner  could  account  for  input  disturbances  and 

“2 

model  uncertainties  which  do  not  depend  upon  0 (t) . The 
-2 

factor  2 0 (t)  would  reflect  operator  adaptability 

throughout  the  trajectory  as  the  operator  perceives  a 
changing  situation, with  permitting  adjustment  of  this 
capability. 

7.3.6  Summary 

In  this  chapter,  the  identif iability  of  the  parameters 
of  the  optimal  control  model  was  analyzed  using  the  theory 


developed  in  the  earlier  chapters.  It  was  found  the  param- 
eters of  the  model  are  identifiable  when  operating  in  a 
tracking  task.  This  is  a very  important  result  from  the 
standpoint  of  the  viability  of  the  model  structure f since 
model  validation  is  dependent  on  the  parameters  being  iden- 
tifiable under  some  combination  of  input  and  measured  out- 
put data.  This  combination  of  input  signal  and  observation 
process  that  yields  identif lability  is  certainly  not  unique. 
Under  certain  observation  processes # some  of  the  parameters 
may  not  be  identifiable  (Ref  40),  but  this  should  not  be 
used  as  a basis  for  claiming  unidentifiability  unless  it  is 
the  only  available  observation  process.  The  input  signal 
(angular  acceleration)  and  observation  process  (ensemble 
averages  and  standard  deviations  of  tracking  errors  for 
many  trials)  was  chosen  to  conform  to  available  simulator 
data.  It  was  fortuitous  that  such  a combination  does  permit 
identif lability . Although  this  approach  differs  from  the 
more  usual  approach  of  basing  identification  on  time  histo- 
ries of  individual  trials,  it  has  the  advantage  of  matching 
model  response  to  the  average  response  of  the  actual  system 
over  many  trials  and  many  operators.  Also,  one  obtains,  as 
a derived  measured  quantity,  the  standard  deviation  of  this 
response,  which  provides  additional  information  for  the 
identification  process.  The  concept  of  using  a deterministic 
trajectory  as  the  basis  for  the  input  signal  (target  angular 
acceleration)  allows  repeatability  of  the  experiment  over 
many  trials  and  forms  a basis  for  predicting  a nonzero  mean 


response.  As  will  be  discussed  below,  the  deterministic 
input  also  would  permit  a systematic  investigation  of  out- 
put sensitivities  to  various  inputs  (trajectories) . 

For  the  computations  that  were  performed,  the  mean 
tracking  error  and  the  tracking  error  standard  deviation 
were  used  as  the  measured  data.  The  results  of  Section  7.2 
indicate  that  this  observation  process  should  provide  suf- 
ficient data  to  permit  identification  of  model  parameters. 
In  the  computations  it  was  decided  to  adjust  each  of  the 

model  parameters  at  each  iteration  using  the  gradient 

Ic 

technique  for  obtaining  the  parameter  change,  A4)  . This 
resulted  in  the  cost  functional  being  reduced  at  each 
iteration,  thereby  verifying  the  technique;  however,  fur- 
ther research  could  test  other  computation  techniques  with 
various  combinations  of  inputs  and  measured  data.  Another 
important  factor  is  the  assessment  of  sensitivity  of  output 
data  to  the  parameter  adjustments  for  various  inputs.  A 
measure  of  output  sensitivity  is  obtained  from 


‘3yi(t^)]2 


where  is  the  i^^ 

y(tj)  and  4.^^  is  the 


component  of  the  measured  vector 
component  of  the  parameter  vector  4>- 


This  sensitivity  function  could  be  evaluated  for  each  param- 


eter under  various  inputs  and  observations.  The  matrix 
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! 

3y(t.)‘ 

3(^ 

yields  some  information  concerning  the  sensitivity  of 
measured  data  to  changes  in  ({>  for  the  given  input  signal. 

In  the  computations  reported  previously,  we  found  that  the 
determinant  of  M((|>  ) was  largest  for  the  trajectory  2 cal- 
culations than  for  trajectory  1,  and  the  smallest  value  was 
obtained  for  the  elevation  axis  of  trajectory  1.  This  can 
be  attributed  to  the  relatively  small  angular  acceleration 
exhibited  by  the  target  in  elevation  for  this  trajectory. 
Thus,  the  output  data  apparently  is  less  sensitive  to 
changes  of  ({>  for  the  elevation  axis  of  this  trajectory  than 
for  the  other  three  cases . More  research  is  warranted  on 
the  subject  of  input  signal  sensitivity. 

The  computational  results  indicate  that  the  optimal 
control  model  does  adequately  predict  human  response  in  a 
tracking  task.  The  mean  response  and  standard  deviation  of 
the  response  of  the  model  predictions  were  in  general 
agreement  with  those  of'^he  simulator  for  each  case  pre- 
sented. As  mentioned  in  the  previous  discussion  of  the 
results,  further  refinements  to  the  model  are  warranted. 

Such  refinements  should  consider  the  possibility  of  cross- 
coupling between  axes  and  methods  for  modeling  such  effects. 
Also  the  effect  of  the  changing  shape  of  the  finite  image  on 
the  display  needs  to  be  researched.  In  Section  6.6,  a 
method  was  presented  for  modeling  the  indifference  threshold 
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to  account  for  finite  image  size;  however,  this  does  not 
account  for  the  possibility  of  the  changing  image  affecting 
the  mean  response.  That  is,  for  a given  trajectory,  the 
image  presented  is  a deterministic  quantity  as  a function 
of  time  which  may  influence  the  operators  in  positioning 
their  sights. 

In  the  next  chapter  we  will  summarize  the  results  of 
this  research  and  recommend  topics  for  further  research. 


Chapter  VIII 


SUMMARY  AND  RECOMMENDATIONS  FOR  FUTURE  WORK 

In  this  chapter,  the  major  objectives  and  contributions 
of  this  research  will  be  reviewed,  and  then  some  potential 
areas  of  extension  emd  future  research  will  be  discussed. 

8.1  Objectives  and  Contributions 

The  first  objective  was  to  develop  sufficient  con- 
ditions for  establishing  observability  and  identifiability 
of  nonlinear  dynamical  systems.  We  considered  systems 
described  by 

x(t)  = f(t,  x(t),  u(t))  x(t  ) = X (8.1.1) 

o o 

for  all  t G C'^Qf  and  an  associated  observation  process 

of  the  foiim 

y(t  ) = h(t  , g(t  , X , u) ) i = 1,  2,  • • • k 

X X X O 

(8.1.2) 

where  x(t)  is  an  n vector,  y(t^)  is  an  m vector  with 
ti  e [t^,  tf]  for  all  i,  u is  an  r vector  valued  function 
of  t,  and  g(t,  x^,  u)  is  a solution  of  Eq  (8.1.1)  for  all 
t e [t^,  tj].  The  observability  problem  is  concerned  with 
determining  conditions  under  which  knowledge  of  the  input 
and  observed  output  data  uniquely  determine  the  state  of 
the  system.  The  approach  taken  was  to  view  the  observation 
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sequence  {y(t^)}  as  an  ink-vector 


y “ Cy’’(ti)y’’(t2)  • • • y'^(tj,)] 


(8.1.3) 


Then  the  nonlinear  system  is  observable  under  the  obser- 
vation process  Y with  the  input  u if  there  is  a one-to-one 
correspondence  between  the  initial  conditions  and  the 
observation  vector  Y(Xq).  The  approach  allowed  us  to  use 
the  inverse  function  theorem  as  a basis  for  stating  suffi- 
cient conditions  for  local  observability  under  discrete 
observation  processes  in  terms  of  the  rank  of  the  Jacobian 
matrix  3Y(Xq)/3Xq  (see  Theorem  2.2).  An  alternate  condition 
based  on  Theorem  2.2,  but  using  the  properties  of  the  Gram 
matrix,  requires  the  nonsingularity  of 

k r3y(t.)-|T  r 9y(t.)-| 

M«^)  = I 


for  local  observability  (see  Theorem  2.3).  Application  of 

these  theorems  require  a minimum  of  k = n/m  observation 

times  and,  in  general,  an  a priori  knovrledge  of  the  solution 

so  that  an  expression  can  be  obtained  for  y(t.)  = h(t.,x  ,u)) 

1 1 o 

An  alternative  test  which  avoids  these  problems  is  derived 
based  on  an  nm-vector  F(t,  x)  where 

P(t,  X)  = [F^(t,  x)F^(t,  X)  • • • X)]  (8.1.4) 


and  F^(t,  x)  given  by  the  recursion  relation  Eq  (2.1.14). 
Using  the  previous  results  in  conjunction  with  the  recursion 
relation,  sufficient  conditions  for  local  observability  are 


I 


Ir 


stated  in  terms  of  the  rank  of  3F(t,  x)/9x  evaluated  at 
t*  e {t^}  (see  Theorem  2.6).  These  results  are  extended  to 
systems  with  continuous  observation  processes  in  Theorems 
2.4,  2.5,  and  2.7. 

Sufficient  conditions  for  global  observability  are 
obtained  for  an  open,  path-connected  set  by  adding  to  the 
conditions  for  local  observability  the  requirement  of  norm- 
coerciveness  or  continuous  differentiability  (see  Theorem 
2.11).  Alternative  sufficient  conditions  for  global  observ- 
ability are  established  using  the  properties  of  strictly 
monotone  and  convex  functions  (see  Theorems  2.14  and  2.16). 
These  results  are  combined  with  the  recursion  relation  of 
Eq  (2.1.14)  to  give  a very  useful  observability  test  (see 
Theorem  2.18).  Again  these  results  are  extended  to  con- 
tinuous time  measurements  in  Theorems  2.17  and  2.19. 

The  results  sximmarized  above  provide  new  and  useful 
conditions  for  establishing  observability  of  nonlinear 
dynaunical  systems.  As  shown  in  Section  2.4,  the  theory 
which  is  developed  can  also  have  useful  applications  to 
linear  systems. 

In  Chapter  III,  the  problem  of  identifying  parameters 
of  nonlinear  dynamical  systems  is  addressed.  Here  we  con- 
sidered systems  described  by 

x(t)  = f(t,  x(t),  u(t),  ({>)  x(t  ) = X (8.1,5) 

o o 

for  all  t e where  is  a vector  of  s parameters. 

The  parameter  vector  of  the  nonlinear  system  is 
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identifiable  under  the  observation  process  with  input  u and 
initial  conditions  if  there  is  a one-to-one  correspond- 
ence between  the  parameter  vector  (|)  and  the  observation 
vector  Y(4i).  By  viewing  parameter  identifiability  as  a 
special  case  of  observability,  sufficient  conditions  for 
identifiability  are  established  based  on  the  development  in 
Chapter  II.  Theorems  3.1,  3.2,  and  3.3  state  sufficient 
conditions  for  local  identifiability  for  discrete  and  con- 
tinuous observations  based  on  the  rank  of  3Y/3<j)  and  the 

k T 

nonsingularity  of  |_^[3y  (t^^) /3(})]  [3y  (t^) /3(J)]  . The  recursion 

relation  of  Eq  (2.1.14)  is  used  to  develop  alternative 
sufficient  conditions  for  local  identifiability  (see 
Theorem  3.4).  Again,  using  the  development  in  Chapter  II, 
sufficient  conditions  for  global  identifiability  on  an  open, 
path-connected  set  are  established  by  adding  to  the  require- 
ments for  local  identifiability  continuous  differentiability 
or  norm-coerciveness.  Also  global  identifiability  is  estab- 
lished based  on  strict  monotinicity  relative  to  4)  of  the 
observations  (see  Theorems  3.6,  3.7,  and  3.8).  Stochastic 
identifiability  is  investigated  for  the  case  of  additive 
zero-mean  white  noise  on  the  observation.  It  was  found 
that  deterministic  identifiability  was  both  necessary  and 
sufficient  for  stochastic  identifiability  for  the  additive 
noise  case  when  independent  repetitions  of  the  observation 
process  can  be  obtained  (see  Section  3.5). 


I 
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The  theory  developed  in  Chapter  III  provides  new  tests 
for  establishing  identif iability  of  nonlinear  dynamical 
systems.  These  results  are  considered  significant  in  view 
of  the  fundamental  nature  of  the  identif iability  question 
in  conjunction  with  system  modeling  and  because  in  many 
situations  the  state  equations  are  not  linear  with  respect 
to  the  parameter  vector  <{i.  As  exemplified  by  Section  3.4, 
these  developments  have  useful  applications  to  linear 
systems,  as  well. 

The  second  major  objective  of  this  research  was  to  use 
the  theory  developed  for  identifying  nonlinear  dynam.ical 
systems  to  investigate  the  identif iability  of  the  optimal 
control  model  for  human  operators.  The  model  equations 
were  developed  to  conform  to  a simulator  where  the  measured 
data  consists  of  ensemble  averages  and  variances  of  many 
trials  of  tracking  data  using  the  simulator.  We  found  that 
each  of  the  parameters  are  globally  identifiable  on  (0,  ®) 
from  observations  of  the  mean  and  standard  deviation  of  the 


tracking  errors  or  of  the  mean  and  standard  deviation  of 
the  manual  control  inputs  to  the  tracking  system  (see 
Section  7.2).  This  is  a significant  result  because  the 
fundam.ental  characteristic  of  parameter  identifiability  had 
not  been  established  for  this  model,  even  though  it  had 


been  successful  in  simulating  human  response.  Further 
significance  is  attached  to  this  result  since  it  exemplifies 
an  application  of  the  theory  to  a nontrivial  nonlinear 
dynamical  identification  problem.  An  additional  result 

217 


O’  r ' 


I 


concerning  the  local  identification  of  the  model  parameters 
was  obtained  from  the  computational  results,  where  we  found 
that  for  each  case  the  matrix 


k 

M(<|))  = I 

i=l 


3y(ti)- 

T 

'ay(t.)' 

3(|) 

3(j) 

was  nonsingular  (see  Section  7.1).  From  Theorem  3.2,  this 
implies  the  parameter  vector  (J)  is  locally  identifiable. 

Computations  were  performed  using  the  gradient  tech- 
nique to  find  an  estimate  of  the  model  parameter  vector  (}> 
which  minimizes  the  cost  functional 


k T 

J = I [z(t  ) - y(t  )]  [z(t  ) - y(t  )] 

i=l  1 1 1 1 


where  the  measured  data  was  the  azimuth  and  elevation  mean 
tracking  errors  and  standard  deviations  (see  Section  7.3). 
Each  axis  was  treated  as  an  independent  case  since  the 
system  had  two  people  in  the  tracking  system;  one  operated 
the  azimuth  control  and  one  the  elevation  control.  Starting 
with  nominal  values  for  the  parameters,  the  computational 
technique  was  successful  in  reducing  the  cost  functional 
with  each  iteration  in  adjusting  the  param.eter  values. 

These  results  are  considered  significant  since  they  estab- 
lish the  feasibility  of  using  a computational  technique  in 
conjunction  with  a cost  functional  for  estimating  the  model 
parameters,  as  opposed  to  previous  work  where  model  param- 
eters were  modified  in  a heuristic  fashion  to  fit  measured 
data. 
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A third  objective  was  to  provide  further  insight  into 
the  modeling  of  human  response  by  the  optimal  control  model. 
Although  the  modeling  approach  was  basically  that  used  by 
Kleinman  (Refs  20^23)  , there  were  several  unique  aspects. 
System  dynamics,  sensor  displays,  and  control  mechanisms 
differ  from  previous  studies;  particularly  the  use  of  two 
person  controls.  In  addition  to  these  unique  aspects,  some 
significant  aspects  of  the  model  which  were  illuminated  by 
this  research  include; 

(1)  An  apparent  operator  indifference  threshold  due 
to  the  finite  image  size  of  the  target  on  the  monitor  which 
increases  the  magnitude  of  the  observation  noise  covariance 
when  tracking  within  the  threshold  (see  Section  6.6).  An 
additional  effect  of  the  changing  image  shape  on  the  miean 
tracking  response  is  suggested  by  comparing  model  and  simu- 
lator responses  (see  Section  7.3). 

(2)  The  results  support  the  model  configuration  which 
treats  target  angular  acceleration  as  a plant  disturbance 
(see  Section  5.4  and  7.3).  In  addition,  the  results  sug- 
gest some  cross  coupling  between  axes  which  could  be 
modeled  by  adding  an  equivalent  observation  noise  in  the 
azimuth  axis  that  is  proportional  to  the  variance  of  the 
observed  quantities  in  the  elevation  axis,  and  vice  versa 
(see  Section  7.3). 

The  above  findings,  together  with  the  general  agreement 
of  the  model  tracking  response  with  the  simulator  tracking 
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There  are  several  aspects  of  the  nonlinear  observ- 
ability and  identification  problem  which  warrant  additional 


effort.  An  area  mentioned,  but  not  explored  in  detail  in 
this  research,  is  the  effect  of  the  input  function  on  the 
observability  and  identification  problem.  There  are  two 
aspects  of  this  phase  of  the  problem  which  could  be 
researched  more  thoroughly.  One  is  the  effect  of  the  input 
function  on  the  binary  question  of  whether  a system  is 
observable  or  identifiable.  The  second  area  concerns  the 
effect  of  the  input  function  on  the  sensitivity  of  output 
data  to  changes  in  initial  conditions  or  system  parameters. 
As  mentioned  in  Chapter  VII,  a measure  of  output  sensitivity 
is  obtained  by  evaluating 


m k 

= I 1 

i=l  j=l 


3yj(tj) 


1 2 


3<J>( 


where  y^(tj)  is  the  i^^  component  of  the  measured  vector 
y(tj)  and  is  the  component  of  the  parameter  vector  (J>. 

The  effect  of  the  input  signal  on  this  sensitivity  function 
needs  to  be  researched  for  nonlinear  systems. 

In  Chapter  III,  a special  case  of  stochastic  identifi- 
ability  was  investigated  for  nonlinear  dynamical  systems 
with  additive  zero-mean  white  noise  on  the  observations. 
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There  is  need  for  further  effort  on  this  difficult  problem 
of  stochastic  identifiability  of  nonlinear  systems. 
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In  Chapters  II  and  111,  the  emphasis  was  placed  on 
establishing  sufficient  conditions  for  observability  and 
identifiability  of  nonlinear  systems.  A better  under- 
standing of  the  problem  would  result  from  a formulation  of 
both  the  necessary  and  sufficient  conditions.  Also,  some 
additional  effort  would  be  useful  in  relating  the  classical 
theory  of  nonlinear  differential  equations  to  the  problem 
of  observability  and  identifiability.  For  example,  how  do 
the  requirements  for  observability  affect  the  properties  of 
the  solution  to  the  system  differential  equation? 

There  are  areas  requiring  further  effort  related  to 
the  identification  of  the  optimal  control  model  parameters. 
Although  the  identifiability  of  the  optimal  control  model 
parameters  was  established  by  the  research  reported  in  this 
document  under  a given  combination  of  input/output  data, 
additional  research  needs  to  be  done.  For  example,  the 
model  was  placed  in  the  context  of  a tracking  task  in  this 
research  to  conform  with  an  available  simulator.  The  iden- 
tifiability of  the  model  parameters  in  other  applications, 
such  as  pilot  modeling,  and  under  other  input/output  com- 
binations needs  to  be  investigated. 

The  parameter  sensitivity  question  needs  to  be  addressed 
in  conjunction  with  identifying  the  parameters  of  the  optimal 
control  model.  This  is  an  important  aspect  of  the  problem  of 
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determining  accurate  estimates  of  the  parameter  values, 
once  identifiability  has  been  established.  As  mentioned 
previously,  although  the  gradient  technique  used  in  this 
research  appeared  to  be  working  satisfactorily,  additional 
effort  is  warranted  to  investigate  the  effectiveness  of 
other  computational  techniques.  Additional  computational 
work  needs  to  be  done  in  conjunction  with  possible  model 
refinements  to  estimate  the  model  parameter  values  more 
confidently. 

Additional  analysis  is  required  to  assess  the  need  for 
further  model  refinements.  Further  effort  is  needed  with 
respect  to  cross  coupling  between  axes  and  target  image 
effects  on  operator  threshold  and  biasing.  Another  area 
needing  investigation  is  the  operator  cost  functional;  for 
example,  during  the  acquisition  phase,  where  tracking  errors 
are  large,  the  mean  squared  error  cost  functional  probably 
does  not  apply. 

This  research  establishes  the  basis  for  addressing  the 
areas  outlined  above.  Clearly,  there  is  a need  for  consid- 
erably more  research  concerning  the  properties  of  nonlinear 
dynamical  systems.  Although  the  nonlinear  aspect  can  slow 
the  progress  of  the  researcher,  the  reward  can  be  great 
because  of  the  extreme  importance  of  this  class  of  systems. 
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Appendix  A 

PRECHET  AND  GATEAUX  DERIVATIVES 

Some  basic  properties  of  Frechet  and  Gateaux  deriv- 
atives of  n-dimensional  functions  will  be  reviewed  in  this 
section  (see  Ref  39:59). 

Definition:  A mapping  F : D C r”  -*•  R™  is  Gateaux-differen- 

tiable  at  an  interior  point  x of  D if  there  exists  a linear 
operator  A e L(r”,  r”*)  such  that,  for  any  h e r” 

lim  (1/t)  II  F(x  + th)  - F(x)  - t A(h)  M = 0 
t-*-0 

If  F is  Gateaux  differentiable,  the  linear  operator  A is 
unique,  is  denoted  by  F' (x) , and  is  called  the  G-derivative 
at  X. 

Definition:  A mapping  F : D C r”  R"*  is  Frechet-dif feren- 

tieUple  at  X G int(D)  if  there  is  an  A e L(r”,  r"*)  such  that 

1 

lim  (1/llhll)  II  F(x  + h)  - F(x)  - A(h)  II  = 0 
h-^0 

If  F is  Frechet  differentiable,  the  linear  operator  A is 
unique,  is  denoted  by  F' (x) , and  is  called  the  F-derivative 
at  X. 

n m 

Definition:  The  F-derivative  of  F : D C R -*■  R at  x^  e D 

is  strong  if,  given  any  e > 0,  there  is  a 6 > 0 so  that  the 
set  S(Xq,  6)  = {y  e r”  3 II  y - | I < 6}  is  a subset  of  D 


JSUJUWUJJlvJJiy 


and 

I I P(y)  - P(x)  - F*  (x^)  (y  - X)  I I < e 1 I y - X I I 
for  all  X,  y e S(x^,  6). 

A.l 

Suppose  that  F has  an  P-derivative  at  each  point  of  an  open 
neighborhood  of  x^  e D.  Then  P'  is  strong  at  x^  if  and  only 
if  F*  is  continuous  at  x. 

A.  2 

P can  have  a strong  derivative  at  x^  even  though  F is  not 
^^^f®J^®^tiable  at  all  points  of  an  open  neighborhood  of  x^, 

A.  3 

If  P : D C r"  -*•  R™  is  F~differentiable  at  x,  then  F is  con- 
tinuous at  X.  That  is,  there  is  a 6 > 0 and  a C > 0 such 
that 

II  F(x  + h)  - F(x)  II  < C I Ihl I 
Whenever  | |h| | ^6. 

A. 4 ^ 

P is  G-differentiable  at  x if  it  is  P-differentiable  at  x. 

A. 5 

P' (x)  is  continuous  at  x if  and  only  if  all  the  partial 
derivatives  3fj^/3xj  are  continuous  at  x. 


r 


A. 6 


If  F has  a G-derivative  at  each  point  of  an  open  neighbor- 
hood of  X and  F' (x)  is  continuous  at  x,  then  F is  F-differ- 
entiable  in  an  open  neighborhood  of  x and  the  F-derivative  is 
continuous . 


A.  7 


If  F has  a continuous  G-derivative  on  the  open  set  D^C  D,  we 


will  say  F is  continuously  differentiable  on  Dq. 


i-i 
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Appendix  B 

CONDITIONS  FOR  UNIQUE  SOLUTIONS  OF  DIFFERENTIAL  SYSTEMS 

In  this  appendix,  we  will  briefly  review  some  con- 
ditions which  will  assure  that  the  system 

x(t)  = f(t,  x(t),  u(t))  = Xq 

for  all  t e t^]  has  a unique  solution  (denoted  by 

g(t,  X , u) ) on  [t^,  t^]  (see  Ref  47). 

Definition;  A mapping  F : D C r”  -►  R™  is  Holder-continuous 
on  Dq  C D if  there  exist  constants  C ^ 0 and  p e (0,  l]  so 
that  for  all  x,  y e D, 

I I F(y)  - F(x)  M<C  II  y-x  II  P. 

If  p = 1,  then  F is  Lipschitz-continuous  on  Dq. 

B.l 

Given  the  system 

X(t)  = f(t,  X(t),  U(t))  x(tQ)  = Xq  (B.l) 

Let  £(•,•,•)  be  piecewise  continuous  in  the  first  argument, 
Lipschitz-continuous  in  the  second  argument,  and  continuous 
in  the  third  argument  for  all  t e j^t^,  t^]].  Then  given  any 
tjj,  Xq  and  piecewise  continuous  u(»),  the  system  differ- 
ential Eq  (B.l)  has  one  and  only  one  solution  g(t,  Xq,  u) 
on  Ct^»  ^f3*  Thus  for  any  t-^  e [t^r  ^f3 » solution 

g(ti,  Xq,  u)  is  uniquely  determined  by  t^,  x^,  and 

"Ct„,  tj- 


226 


Appendix  C 


EXISTENCE  OF  SOLUTION  TO  RICCATI  EQUATION 


We  want  to  show  that  our  formulation  of  the  linear 
regulator  problem  results  in  the  existence  of  a solution  to 
the  steady  state  Riccati  equation. 

Consider  the  linear  equation 


x(t)  = A x(t)  + B u(t) 


(C.l) 


with  observations 


Z(t)  = D x(t) 


(C.2) 


and  a cost  functional 


J = / [z^(t)Q  Z(t)  + g u^(t)]dt 

o 


(C.3) 


with  Q positive  definite  and  g > 0. 

The  associated  Riccati  equation  for  the  linear  regulator  is 
- P(t)  = P(t)A  + A*^P(t)  - P(t)B  g”^  B*^p(t)  + d'^’q  D (C.4) 


Kwakernaak  and  Sivan  (Ref  10)  have  stated  conditions  for  the 
Riccati  equation  to  have  a steady  state  solution.  As 
tj  -*•  the  solution  of  the  Riccati  equation  approaches  a 
constant  steady  state  value  if  and  only  if  the  system  as 
defined  by  Eqs  (C.l)  and  (C.2)  possesses  no  poles  that  are 
at  the  same  time  unstable,  uncontrollable,  and  reconstruct- 
able 


In  the  context  of  the  optimal  control  model  in 


Chapter  V,  we  have 


r^n' 
1 

0 1 
1 

0 * 

A = 

^d  1 

®c 

1 

Lo  , 

0 1 

0 . 

’ 0 ' 

f X 1 

B = 

0 

n X 1 

1 , 

1x1 

D = 

! 1 “J 

1x1 


Ajj  = n X A 

A^  = n X n 

B = n X 1 

c 


D = (il  + n + l)x(Jl  + n + 1) 
= (Jl  + n + 1)  X n 


Q,  = n X n 
d 


with  the  above  forms  for  the  matrices  A,  B,  D,  and  Q it  is 
clear  that  the  first  i states  (disturbance  states)  are 
unstable  and  uncontrollable.  However,  to  be  reconstructable, 
the  present  state  x (t)  must  be  determined  from  knowledge  of 
past  outputs  {z  (a)  = D x(a)  : o £ t).  From  the  form  of  D, 

It  is  clear  that  the  first  i states  are  not  reconstructable. 
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rv^  V . 


fti*' 


Although  the  states  of  the  dynamics  (states  A + 1 through 
A + n)  may  be  reconstructable » they  are  also  controllable. 
Thus  our  formulation  of  the  linear  regulator  problem  leads 
to  the  existence  of  a solution  to  the  steady  state  Riccati 
equation . 
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Appendix  D 


EFFECT  OF  HOLDING  NEURO-MUSCULAR  DELAY,  T , 
CONSTANT  ON  SOLUTION  OF  RICCATI  EQUATION 


In  Section  5.3,  it  was  discussed  that  there  is  physical 
motivation  to  hold  the  neuro-muscular  delay,  t^,  constant. 

We  want  to  show  that  for  the  class  of  problems  under  con- 
sideration, specifying  a value  of  completely  determines 
the  optimal  control  gains,  X. 

Recall  from  Eq  (5.3.9)  that  the  optimal  gains  are 
obtained  by 


g [^13  I **23  1 ^33] 


Let  A 


13 


P Po^ 

i-j  23 

and  = p — 
P33  23  P33 


so  that 


^13  I ^23  I ^ 


] 


From  Eq  (5.30)  we  know ^ 


* p 


33 


Thus  Eq  (D.2)  becomes 


^ '‘13 1 ‘23 : ^ 


] 


(D.l) 


(D.2) 


(D.3) 


(D.4) 
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We  will  show  below  that  for  the  cost  functionals  with 
weighting  matrix 


i.  and  S,  can  be  obtained  independent  of  the  value  of  q. 
13  23 

That  is,  the  value  of  completely  specifies  X. 

With  constant,  the  partitioned  Riccati  equations 
(5.3.5),  (5.3.6),  and  (5.3.7)  become 


T Qd 

^d  ■'^n^d  *"22 


+ A 


2,  =0 

23 


With  this  form,  the  factor  q appears  in  only  one  of  the 
above  system  of  [n(n+l)/2]  + n + 1 equations.  If  we  ignore 
the  equation  involving  q,  then  the  system  represents 


I 

\i 


[n(n+l)/2]  + n equations  with  the  same  number  of  unknowns 
(i.e.,  n<n+l)/2  components  of  I22  components  of  ^23^* 

Thus  we  can  solve  for  the  values  of  I22  ^23 

of  the  value  of  q.  The  equation  where  q appears  determines 
the  value  of  P33  terms  of  ^22'  ^23'  snd  q. 

Now  define  ^11  “ ^11/^33'  *"12  ~ ^12^^33' 

= Pq^3/P33*  Then  equations  (5.3.2),  (5.3.3),  and  (5.3.4) 

become 


T iiiA^  + T ^ A,  + T A S,_-  + T S',  ~ ~ 0 

11^  n 12  D n n 11  n b 12  13  13 


(D.8) 


tJIA+ta'*^^  +t  a * 0 

n 12  d n n 12  n b 22  13  23 


(D.9) 


+ A - T Jl  =0  (D.IO) 

12  c n 13  b 22  n 13 


With  i.  - and  A known,  the  eibove  system  represents 

AM  A ^ 

CA(A+1)/2  + (Axn)  + a]  equations  and  the  same  number  of 
unknowns . 

From  the  ed>ove,  it  is  apparent  that  for  a specified 

value  of  T_  we  can  solve  the  partitioned  Riccati  equation 
n ■%_:> 

for  A33  and  A23  independent  of  the  value  of  q.  Then 


i - T„  I *23  1 • 


(D.ll) 


Appendix  E 

SIMULATOR  RICCATI  EQUATION  SOLUTION 

Cy2  Y3  Y4^ 

ys 

^*12  ^13  *14^ 

*22  *23  *24' 

*32  *33  *34 

.*42  *43  *44. 

Eqs  (6.4.5)  through  (6.4.10)  become: 


Let  P 


13 


23 


33 


T A 
J 

T ^ 
A 

A 


12 


22 


0 


*12  *13  *14^ 


cos0g 


+ [0  0 -1]  [ Y^  ] 


“ t ^1  ^5  = ° 


(E.3) 


*22  *23  *24 


*32  *33  *34 


*42  *43  *44J 


0 1 
-a,  -e 


1 "1 
10  0 


0 -a. 


1 -6, 


*22  *23  *24 


*32  *33  *34 


*42  *43  *44 


o 

o 

o 

^2 

0 0 O' 

+ 

0 0 0 

- 1 ^^2  ^3  ^4^ 

^3 

= 

0 0 0 

.0  0 q 

.^4 

o' 

O 

O 

(E.4) 


*22 

*23 

*24 

0 

0 

-a 

1 

1 

to 

j 

*32 

*33 

*34 

cos0„ 
2 E 

+ 

1 

-^1 

0 

^3 

*42 

*43 

*44. 

0 

0 

0 

0 

.y4. 

1 

g 


^2 

O' 

y3 

ys  = 

0 

^4- 

.0 

(E.5) 


[^2  ^3  ^4] 


Kz  cosGg 


= 0 


(E.6) 
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Eqs  (E.4),  (E.5),  and  (E.6)  can  be  solved  for  the  control 


lable  state  gains.  Expanding  Eqs (E.4) , (E.5),  and  (E.6), 
one  can  obtain  the  following  set  of  equations; 

2(-  X32  + X42)  = I ^12  (E. 


^2  - ^1^32  = g ^4^3 

(E.l 

^^^23  “ ^1*33^  g ^3 

(E.l 

Y2X33  + Y2  - 6^73  ' ^ Ya  » 0 

(E.ll 

(E.i; 

^2^3 

Xj2  ®x*23  “1*33  ^43  " g 

(E.i: 

^2*23  - “l^3  ^4  - 9 " “ 

(E.i: 

y5Y4 

^2''43  ■ g " ° 

(E.l^ 

1 2 
^^4=^5 

('E.l! 

where  y,  “ E cos©  . 
z 2 E 

The  above  nine  equations  have  nine  unknowns  which  can  be 
solved  as  outlined  below. 


From  these  one  can  show  that 


V4 


Ys  = N2  g ^2  y's 


(E.16) 

(E.17) 


Ml  - 4 '‘3^1  ] 


(E.18) 


where 


2 ^5^3 

El  = A^y3  + AgYi  + AgYs  + — 


^3^2 


2 ^ + 

3i9 


-Y^y 


2 4 


^3  ” ^2^5 


(E.19) 


_ -i  = 0 (E.20) 

3, 


where 


A ==^1- 
^2  a^g 


^2 


A3  = 2aj^B^g 


_Il_ 

N ' 2Pig 


"A*Y- 

n ^ ^ I P 


*6  ' - ^ 


The  above  equations  can  be  solved  for  Y2'  ^^3'  ^5 

follows.  Note  that  Y4  is  determined  from  Eq  (E.16).  First 

solve 


1-4  = 0 


(E.21) 


for  maximum  value  of  y^' 
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Note  that  Eq  (E.21)  is  only  a ftinction  of  when  Eq  (E.17) 
and  Eq  (E.19)  are  combined.  This  is  solved  using  the 
Newton  iteration  method. 


®1  " ^^3  “ ° 


* *5=^3  * *6^5 

. ^5^4 

^ g ■ 

1 

■sj 

(E.22) 

*6®  ^2 
+ A_  + 

^ li 

(E.25) 

5 ^5 

^5  ^3 

g 

Using  the  Newton  method 


(E.24) 


This  gives  the  maximum  value  for  y^r  denoted  as  then 

use  y^*  - e as  the  initial  value  for  solving  Eq  (E.20). 
That  is 


s'5^2  ^ r:2£i . !ii 

961  * [SiV  e^J  ''4 


P4  = *3X2  + 


A2Y5 


(£.25) 


^ ^ Yj  dYj  p-YjY^  Oj.| 

= dy3  = 2Ajy3  dYj  + 5^  d^T  ^ 


r ^2  1 

-A  + -4-  — ^ (E.26) 

L 2 g6jJ  yg 


This  gives  values  for  Yy  Y^» 
following  logic  can  be  employed. 
Prom  Eq  (E.2)  one  obtains 


Prom  Eq  (E.4)  one  obtains 


Prom  Eq  (E.14) 


Combining  Eqs  (E.29),  (E.30),  and  (E.31)  we  get  for  y 


Appendix  F 


GRADIENT  EQUATIONS  WITH  RESPECT  TO  PARAMETER  VECTOR  ((> 


F.l  Gradient  with  Respect  to 


3x_  3x  3A  A T At  fle,  3e. 

_ am  ^ ^ u e u e + e " iSl  + 2m 


A T 3e,  3e. 


A 

3x 

am 

’TT 

n 


X . ?.!!gSL  11^  C e + K C 


A = A - B X 
c a a c 


3A  3X  SB 

§.  - B £ - a X 

^ ' 


3A^  r 0 0 


^a  = 

’■'ll  Lo 


-Vt, 


8X^  ax  3, 

^ -5?5v 


Prom  Eq  (5.3.9) 


^ci  “ “ SP  Pi 


(see  Appendix  E) 


for  1 = 1,  2,  3,  4 


(P.8) 


^5  ' • g y.  “ ? g Y,y, 


(P.9) 


Thus  g = 


"n  ^5 


where 


(P.IO) 


yc  + — T 
5 3g  n 


- y . T ^ 

I 22  ^22  ®9J 


Prom  Eqs  (P.8)  and  (F.i 


*oi  “ yj  ''i 


i = 1,  2,  3,  4 


^ = i-  """i  - -J-  Tv  . sy,! 

^ ^5  ^ ® 5y^J  ’'i 

To  evaluate  ^ , „e  use  the  equations  developed  1. 
Appendix  E for  These  are: 


= i-  ^ 


(Fell) 


(P.12) 


(F.13) 


^^2 


2 y 


+ y + f~^2^5  + “ll 


y. 

^3  ■ Vs  ' iq-  ’ “ <'■•14) 


Taking  the  partial  of  Eq  (F.14)  with  respect  to  g yields 


ay,  _ ^ ^ 

”5?  a g6,  9g  „ _„2r 


r '^2^2 


= Ml  + l2i  !l3 


(F.19) 


Y4  = Ng  q 


3y4  _a_ 


Substituting  into  Eq  (F.18)  and  rearranging 

r,  A V + ^1  . [ M M + / IMi  + 

\}  3^2  Bj^gJ  9g  Bj^g 

_ A Ml  ^II  = ^2  y ^ + ^22^2 

2 J 3g  2 a^gB^  2 B^q  B.g^ 


(F.20) 


From  Eq  (F.16) 


^2^3'^  ^2y4y3  ^ 

“l9Bi2y4  „^g28^  2 

^5'^  Vs  , _a_ 

^'’“1^5  9^»X 


^ l—^]ll  . \K 

’9  2Aj2  \ 20jg^8j  / I 


-4A3E1 


(F.21) 


5X7  ( )(  -■'''3  ^ —\  ] 

3 \ 2 1-AA^E^I\  2 a^g^Bi  ' 


(F.22) 


with 


=2a, s.  2“iei9y4 


Substituting  Eq  (F.26)  into  Eq  (P.21)  allows  evaluation  of 


in  terms  of  known  quantities 


From  Eq  (F.17) 


is  given  by  Eq  (F.26)  and 


The  above  allows  evaluation  of 


K = e G 


(F.28) 


(F.29) 


T -1 

G = V 


(F.30) 


(F.31) 


From  the  Eq  (2.3.10)  for  propagating  we  obtain 


3P,  3A^  3P-  3A  3P,  m 

1 = ^P.  + A —i+P,  —a  + ^ A ^ 

Tf;;  axp  1 a 3t„  1 3t„  8t„  a 


3Pi  T -1  T “1  _ . 

C V CP-P.  C V 
a alia  a 


STa  T ‘^*a 

+ ^ wr  + r W ^ 

3Tn  a a 3t 


(F.32) 


0 “1/T. 


(F.33) 


From  Eq  (5.5.2)  we  get 


3e. 

Im 

3A 

i.  e + A 

3e 

3-^n  " 

3Tj^  Im  f 

(F.34) 


. fc’* 


A-  = -•  G C 

fa  a 


3Af  3A 


a 9G 


3t  a 

n n n 


3e 

3t 


T 3Aa  A,e 

e e ^ de 

O l_ 

0 


] r wjCo) 


3x 


The  above  equations  allows  calculation  of  . as  a 


3t 


n 


of  time  using  Eq  (F.l). 


3P. 


To  calculate  — !L  we  refer  to  Eqs  (5.4.27),  (5. 

3-^n 


(5.4.30),  (5.4.31),  and  (5.4.32). 

T 


3^a  _ ‘’"a 

3"n  “ " 3"n 


3 A,  T A^t  gp 


A '^T 


e ^2^  ® 


3t 


n 


A T 


3A  T A '^T  3P,  A 


2 3t 


n 

T * T, 


n 


9A  ^ A^^t  3A  At  „ 

+ ^ e = + X e ^ P 

5 3t„  3 


A.t  31-3  ’'’4  . 

^ ® 3^n  ^ 3t„  ^ 


with 

3P,  3Aa 


^ 3^n  2 


+ G V 


3P, 

3P 

A,  2 

+ 

^ 9^n 

3Tj 

3gT 

2V„ 

m 

3t 

3 

n 

, 3A 

2 W 


n n 


VG 


0 0 
0 1 
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(F.35) 


1 


(F.36)  I 

function  ] 


4.29), 


4 


(F.37) 


(F.38) 


. »■*  V ) 


3P 


3P 

3F  m 

+ A 3 

c 3t^ 

+ __3  A-  + P, 
3t  f 3 

n 

n 

n 

3K 

+ ^ C- 
n “ 

[P2  - Pi]  + K 

3A, 


3 3t 


3P 

3t 


n n 


3P> 


3A 


3t 


n " "n 


3A  3P  3P  T 

c p + A 1 + i A + P 

3t-.  4 c 3t_  3t^  c ^ 3t 


c 12L_  c P '■ 
— + 3t 


n n 


n 


+ K C_  ^ 


* *''3  c ’■  / + P C ’■ 


3K^ 


+ I7-VK’' 

n 


(F.40) 


n 


3Pc  T 3A  A e 
5 f _ a 

^ ' “'o  ^ 

no  n 


a e a p 


0 1 


LO  V„J 
m 


A^e 


T a 
e de 
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2Vm  f ^ a 


T 3 0 

n 


0 

0 


0 

IJ 


e ® de 


T A e 
+ / e ® 

0 


Wi  0 


0 V. 


m 


r/  t e ° dE  (F.41) 


n 


3P. 


The  above  relations  allow  computing  as  a function  of 


n 


time. 

F. 2 Gradient  with  Respect  to  x 


3x 3x 

"St 


3e. 


TT 


A_t 

'a  ~ 'Im  3t 


am  _ ~"am  ^ e ^ e + 


(F.42) 
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aj-«..garjTiraL:-ay  V : 
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(P.43) 


^ ^ is  ^ 

3x  3x 

&TQ  ss  ft  dm  _L.  V K ^ ^ 

TT  c W a Im 


|2i  = A K 

3t  a 


3eo 

2m  a 
— = ® 


T r r' 

.oJ  ®T 


A^T  A_^t  At  T A.^t 

A^  e ^ P2  e ^ + e ^ ?2  Aa  e a 

+ P3  e ® + A e ^ P,  + e ^ 

•»  a a 3 3t 


AaX  ^ 3P4  ?F3 

® 3t  + 3f-  + ^r 


3^3 

= A 

3P3  T 

■*■  A 

+ A K C fp  - 

P 

3t 

c 

3t 

3t  f 

a aL  2 

1 

^^4 

^**4  T 

+ i A 

T 

3P 

3P 

3t 

— A 

C 

3t 

3t  \ 

^a  3t  ^ 

3t 

A '^T 

T T T T i»  ^ T 

\ ""a  ^3  + P3  G \ e + A^  K V K 


A Tt 

T T a 
+ K V G A,  e ^ 
a 


r T a 
^a  ® 


A ^T 


— — ' — — - »■.-•.  ...-»^.»nrr.  -I.-  g.-  _■  ^ ..■■>. 


F.3  Gradient  with  Respect  to  p 


Recall  that 


y(t)  = C X (t)  + v(t) 
« a 


where  v(t)  is  observation  noise  with  covariance 


■Pl  0 1 T ^1  ° 

V(t)  = C P(t)C  = 

0P2  a [O  V. 


In  the  present  case  the  operator  observes  angle  error  and 
angle  error  rate. 
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F.4  Gradient  with  Respect  to  pj^ 

Recall  that  a motor  noise  v^^^Ct)  is  added  to  the  com- 
manded control  u_  with  a covariance 
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where  Pj^^  is  a constant. 
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F.5  Gradient  with  Respect  to  W 


component  of  the  input  disturbance.  The  equations  for  this 
case  is  the  same  as  with  respect  to  p except 


Appendix  G 


Digital  computer  programs  were  developed  to  propagate 
the  mean  states  and  covariances  given  by  the  coupled  matrix 
differential  Eqs  (5.4.24)  to  (5.4.26),  (5.4.12),  and 
(5.4.28)  to  (5.4.31).  Also  the  gradient  equations  developed 
in  Appendix  F for  calculation  of  3J/3(|>  are  incorporated 
into  the  programs.  These  programs  used  the  computational 
algorithms  and  subroutines  developed  by  Klelnman  for  use  in 
Linear  System  Studies  (Ref  24) . Figure  G.l  shows  a flow 
diagram  for  these  programs. 
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Figure  G.l 


Set  Values  of  Parameters 


Set  Initial  Conditions 
Subroutine  Init 


Calculate  Angular  Position  Velocity 
and  Acceleration  for  time  t 

Subroutine  Angle 


Calculate  Feedback  Gains  for  time  t 
Subroutine  FBGl 


Propagate  Coupled  Matrix  Differential 
Equations  for  Mean  States,  Covariances < 
Gradient  Equations  from  time  t 
to  t + At 

Subroutines  Set,  Step,  F2,  F3, 

Dell,  and  Grady. 
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